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Interest  in  the  numerical  solution  of  acoustic  inverse  scattering  problems  arises  in  a  number 
of  areas.  Examples  include  medical  diagnostics,  non-destructive  industrial  testing,  geophys¬ 
ical  prospecting  for  petroleum  and  minerals,  and  detection  of  earthquakes. 

The  highly  nonlinear  and  oscillatory  nature  of  the  problem  is  one  of  the  major  difficulties 
one  encounters  in  the  construction  of  effective  inversion  algorithms.  Schemes  based  on 
global  or  local  linearization  methods,  or  nonlinear  optimization  techniques,  tend  to  work 
only  when  the  index  of  refraction  is  almost  constant.  They  develop  serious  convergence 
problems  whenever  the  perturbation  of  the  index  of  refraction  increases. 

Limited  successes  in  the  solution  of  the  inverse  problems  have  been  achieved  only  in  one  di¬ 
mensional  cases  (Gelfand-Levitan  and  layer  striping  methods  are  among  the  most  notable). 
These  methods  are  generally  unstable  numerically  since  the  procedures  used  to  calculate 
the  index  of  refraction  are  ill-conditioned. 

We  present  a  method  for  the  solution  of  inverse  problems  for  the  one  dimensional  Helmholtz 
equation.  The  scheme  is  based  on  a  combination  of  the  standard  Riccati  equation  for  the 
impedance  function  with  a  new  trace  formula  for  the  derivative  of  the  index  of  refraction, 
and  can  be  viewed  as  a  frequency  domain  version  of  the  layer-stripping  approach.  The 
principal  advantage  of  the  procedure  is  that  if  the  scatterer  to  be  reconstructed  has  m  >  1 
continuous  derivatives,  the  accuracy  of  the  reconstruction  is  proportional  to  l/am,  where  a 
is  the  highest  frequency  for  which  scattering  data  are  available.  Thus,  a  smooth  scatterer 
is  reconstructed  very  accurately  from  a  limited  amount  of  available  data. 

The  scheme  has  an  asymptotic  cost  0(n2),  where  n  is  the  number  of  features  to  be  recovered 
(as  do  classical  layer-stripping  algorithms),  and  is  stable  with  respect  to  perturbations  of 
the  scattering  data.  The  performance  of  the  algorithm  is  illustrated  by  several  numerical 
examples.  Generalizations  of  this  approach  in  two  dimensions  are  discussed. 
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Chapter  1 
Introduction 


1 . 1  Background 

During  the  last  several  decades,  the  inverse  scattering  problems  for  the  Helmholtz 
equation  have  enjoyed  a  remarkable  degree  of  popularity,  both  in  pure  and  ap¬ 
plied  contexts  (see,  e.g.,  [1],  [2]).  A  number  of  algorithms  has  been  proposed 
for  the  numerical  treatment  of  these  problems,  in  such  environments  as  medical 
diagnostics,  non-destructive  industrial  testing,  anti-submarine  warfare,  oil  explo¬ 
ration,  etc.  In  the  design  of  such  a  scheme,  three  major  problems  have  to  be 
overcome. 

1.  The  problem  is  highly  non-linear,  even  in  its  purely  mathematical  form.  In 
the  one-dimensional  case,  the  problem  can  be  reduced  to  a  linear  one,  but  the 
procedure  is  not  stable  numerically. 

2.  Once  a  mathematically  valid  inversion  scheme  is  constructed,  it  might  or 
might  not  be  stable  numerically.  In  fact,  no  numerically  robust  schemes  seem  to 
exist  at  this  time,  except  in  one  dimension. 

3.  The  cost  of  applying  the  scheme  on  the  computer  tends  to  be  extremely  high, 
except  in  the  one-dimensional  case. 

The  existing  attempts  to  solve  inverse  scattering  problems  for  the  Helmholtz 
equation  can  be  roughly  subdivided  into  four  groups. 

1.  Linearized  inversion  schemes,  attempting  to  approximate  the  inverse  scattering 
problem  by  the  problem  of  inverting  an  appropriately  chosen  linear  operator  (see, 
for  example,  [2]). 

2.  Methods  based  on  the  non-linear  optimization  techniques,  attempting  to  re¬ 
cover  the  parameters  of  the  problem  iteratively,  by  solving  a  sequence  of  forward 
scattering  problems  (see,  for  example,  [3],  [4],  [5]). 
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3.  Gel’fand-Levitan  and  related  techniques,  converting  the  Helmholtz  equation 
into  the  Schrodinger  equation,  the  inverse  problem  for  the  latter  being  reducible 
to  the  solution  of  a  linear  Volterra  integral  equation  (see,  for  example,  [1],  [6]). 

4.  Techniques  based  on  the  so-called  trace  formulae,  connecting  the  high  fre¬ 
quency  behavior  of  the  solutions  of  the  Helmholtz  equation  with  the  local  values 
of  the  parameters  to  be  recovered  (see,  for  example,  [7],  [8],  [9]). 

The  approach  of  this  thesis  falls  into  the  category  4  above,  and  is  different 
from  the  preceding  work  in  the  choice  of  the  trace  formula  (see  Theorem  4.19 
in  Section  4.3).  The  new  trace  formula  leads  to  an  algorithm  with  superior 
convergence  properties  for  smooth  scatterers  (see  Chapter  5  below),  and  the 
resulting  numerical  procedure  is  extremely  stable  and  efficient. 


1.2  Thesis  Organization 

Chapter  2  contains  the  exact  formulation  of  the  problem  to  be  addressed,  to¬ 
gether  with  the  relevant  notation.  In  Chapter  3,  we  summarize  the  background 
facts  to  be  used  in  this  thesis.  Chapter  4  is  devoted  to  the  development  of  the 
mathematical  apparatus  used  to  construct  the  algorithm,  and  in  Chapter  5  the 
scheme  itself  is  presented.  In  Chapter  6  we  present  several  numerical  examples 
demonstrating  the  actual  performance  of  the  procedure.  Finally,  in  Chapter  7 
we  discuss  generalizations  of  the  approach  to  higher  dimensions. 


Chapter  2 

Formulation  of  the  Problem 


2.1  The  Helmholtz  Equation 

Following  the  standard  practice,  we  will  be  considering  the  one  dimensional  scalar 
Helmholtz  equation 


<j>"(x,  k )  +  k2{  1  +  q(x))<f>(x ,  k )  =  0.  (2.1) 

Unless  specified  otherwise,  we  will  be  assuming  that  q  €  cg([0, 1]),  i.e.,  that  q  is 
twice  continuously  differentiable  everywhere,  and  that  q(x)  =  0  for  all  x  £  [0,  lj. 
Defining  the  function  n  :  R  —>  R  by  the  formula 

n(x)  =  y/l+q(x),  (2.2) 

we  will  denote  by  no,  nx  the  minimum  and  maximum  of  n  respectively,  and 
assume  that  0  <  n0  so  that 


n0  <  n(x)  =  yj  1  +  q{x)  <  nx.  (2.3) 

For  any  complex  fc,  we  consider  solutions  of  the  Helmholtz  equation  <£+(2,  k )  and 
(f>-(x,k)  which  have  the  form 

*^+(x,  k)  fiinc+ 1  <f> acat+ (x  1  k^  (2-4) 

<^_(x,/j)  —  <f>inc—(xik)-\-<t>3cat-(x,k)  (2-5) 

with 

<t>inc+{x,  k)  =  e'kx,  (2.6) 

^nc-(x,fc)  =  e"**  (2.7) 
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and  <t>scat+,4>scat-  both  satisfying  the  outgoing  radiation  boundary  conditions 

#e«t(M)  +  **&arf(0,*)  =  0,  (2.8) 

</>'scat(hk)-ik<l>aeat(  l,fc)  =  0.  (2.9) 

Normally,  <j>inc+  and  4>inc-  are  referred  to  as  right-going  and  left-going  incident 
fields  respectively,  and  <pscat+  and  4>scat-  are  called  scattered  fields  corresponding 
to  the  excitations  <&„c+  and  0,„c_.  The  sum  of  an  incident  field  and  its  corre¬ 
sponding  scattered  field  is  called  the  total  field. 

Remark  2.1  Throughout  this  thesis,  given  a  function  f(x,k),  we  will  take  the 
liberty  to  denote  by  f'(x,k),  so  that  the  derivatives  in  the  formulae  (2.15), 
(2.16)  are  with  respect  to  x. 

As  is  well-known,  for  any  complex  k ,  the  scattered  fields  <pscat+(x,k)  and 
4* scat— (x,  k)  satisfy  the  nonhomogeneous  Helmholtz  equations 

4>"scat+{x->k)  +  ^2(1  +  «(*))&«»*+(*>*)  =  -k2q(x)e,kx,  (2.10) 

</>"cat-(xi  k)  +  £2(1  +  q(x))<f>scat-(x,  k )  =  —k2q(x)e~lkx .  (2.11) 

Since  q(x)  —  0  for  all  x  (0, 1),  it  is  easy  to  see  that  for  any  k  £  C  there  exist 

two  complex  numbers  p+(k),  p_(k),  identified  as  the  reflection  coefficients,  such 

that 


4>scat{x,k )  =  p+(k)  ■  e  lkx,  for  all  x  <  0,  (2-12) 

<t>scat{x,k)  =  p-(k)  ■  elkx.  for  all  x  >  1,  (2.13) 

due  to  (2.10),  (2.8)  and  (2.11),  (2.9)  respectively. 


2.2  The  Impedance  Functions 

Denote  by  C+  the  upper  half  of  the  complex  plane  so  that 

C+  =  {k  €  C\Im{k)  >  0}. 

(2.14) 

For  any  k  €  C+,  the  impedance  functions  p+(x,  k),p„(x 
4>+(x,  k ),  <t>-(x,  k),  respectively,  axe  defined  by  the  formulae 

,  k )  associated  with 

P+{  ’  }  ik4>+(x,k)* 

(2.15) 

P-(X'k)  = 

(2.16) 

p-(x,k)  = 
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Remark  2.2  For  x  outside  the  scatterer.  it  is  easy  to  obtain  explicit  expressions 
/orp+,p_  in  terms  of  reflection  coefficients  p+, p_.  Indeed,  combining  (2-4)  with 
(2.12),  (2.5)  with  (2.13),  we  have 


<t>+(x,k)  =  etkx  +  p+(k)e  ,kx,  for  all  x  <  0,  (2-17) 

<f>-(x,k)  =  e~'kx  +  p-(k)elkx,  for  all  x  >  1,  (2.18) 

which  can  be  reformulated  as 

<p+{x,  k)  =  eikx  +  b+{k)e-ikx+Q+W,  for  all  x  <  0,  (2.19) 

d>-{x,k)  =  e~ikx  +  b-(k)eikx+a~(k\  for  all  x  >  1.  (2.20) 


with  a+(k),oc-(k)  real  numbers  and  b+(k)  >  0.  b-(k)  >  0,  for  any  k  €  C . 
Consequently, 


1  —  b\(k)  +  i2b+(k)  sin {kx  —  a+(k)) 
1  +  b\(k)  +  2  b+(k)  cos  (kx  —  a+(k)) 


for  all  x  <  0,  and 


1  —  b2_ (k)  +  i26_ (Ar)  sin(Ax  —  a_(A)) 
1  +  tit(k)  +  2  b-(k)  cos  (kx  —  a-(A)) 


for  all  x  >  1 . 


(2.21) 


(2.22) 


For  any  complex  number  k,  the  boundary  value  problems  for  <f>+,  <t>-  can  be 
reformulated  as  initial  value  problems.  More  specifically,  formulae  (2.4),  (2.5), 
(2.12)  and  (2.13)  imply  that  there  exist  such  complex  constants  a ,  0,  depending 
only  on  k ,  that 


< b+(x,k )  =  a-etkx,  for  all  x  >  1,  (2.23) 

<j>-{x,k)  —  0  ■  e~lkx,  for  all  x  <  0.  (2.24) 

Furthermore,  a,  0  are  nonzero  because,  e.g.,  if  0  =  0,  then  <£_(0,  k)  =  d>'_(0,  k)  = 
0,  according  to  uniqueness  theorem  on  initial  value  problems,  <P-(x,k)  =  0  for 
all  x  €  R,  i.e., 

<t>scat-{x ,  k)  =  —<f>inc-(x,  k)  =  —e~,kx,  (2.25) 

contradicting  to  (2.13).  Clearly,  formulae  (2.23),  (2.24)  can  be  used  as  initial 
conditions  for  equation  (2.1)  to  (uniquely)  determine  the  total  fields  d>+,  <f>_. 

Remark  2.3  While  the  existence  and  uniqueness  of  the  functions  cf>+,4>-  are 
quite  obvious  for  any  complex  k,  the  functions  p+(x,k).p~(x,k)  are  only  well- 
defined  when  Im(k)  >  0,  and  the  proof  of  this  fact  is  somewhat  involved  (see 
lemmas  in  Section  f.l  below). 
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Remark  2.4  It  is  easy  to  see  that  the  impedance  functions  p+,p~  are  indepen¬ 
dent  of  the  nonzero  coefficients  a,  (5  in  (2.23),  (2.24).  Therefore,  for  simplicity, 
the  initial  conditions  (2.23),  (2.24)  are  reformulated  as 

4>+(x,k)  =  eikx,  for  all  x  >  1,  (2.26) 

4>-{x,k)  =  e~ikx.  for  al  x  <  0  (2.27) 

The  functions  as  solutions  of  equation  (2.1)  subject  to  boundary  condi¬ 

tions  (2.26),  (2.27)  differ  from  those  subject  to  boundary  conditions  (2.23), 
(2.24)  by  constants. 

The  classical  inverse  scattering  problem  for  the  equation  (2.1)  is  as  follows: 
Problem  1.  Given  the  impedance  function  p+(0,k)  for  all  k  G  R,  reconstruct 
the  potential  q  for  all  x  €  [0, 1]. 

It  is  well-known  that  this  problem  has  a  unique  solution  (and  in  the  class 
of  functions  q  much  broader  than  Cq([0,  1]),  and  several  constructive  schemes  for 
that  purpose  have  been  proposed,  most  notably  the  Gel’fand-Levitan  and  related 
methods.  However,  in  applications  the  impedance  function  p+( 0,  k )  is  measured 
with  a  finite  accuracy  and  at  a  finite  number  of  (usually  equispaced)  values  of  the 
wavenumber  k.  Therefore,  the  following  problem  is  more  relevant  in  numerical 
applications 

Problem  2.  Suppose  that  the  impedance  function  p+(0,  k)  is  given  at  a  finite 
number  of  frequencies  kj,j  =  1, 2, . . . ,  N  defined  by  the  formulae  kj  =  j  •  h,  with 
h  a  positive  constant.  Suppose  further  that  the  values  p+(0,  k^)  are  given  with 
the  relative  accuracy  e.  Reconstruct  the  potential  q  in  the  interval  [0, 1]  with  the 
error  that  rapidly  decreases  with  increasing  N  and  decreasing  h. 

This  thesis  is  devoted  to  the  construction  of  an  algorithm  for  the  solution  of 
the  Problem  2. 


Observation  2.5  The  value  of  impedance  function  p+  at  x  —  xo ,  x0  <  0  can  be 
obtained  from  <fi+(xo,  k)  in  the  following  manner.  Assuming  that  at  x  <  0,  the 
total  field  <f>+(x,k)  is  given  by  (2.17),  from  which  p+(k)  can  be  obtained 

/<+(*)  =  (#+(*«,  *)  -  eita)  (2.28) 

the  value  of  the  impedance  function  p+  at  x  —  xo  is  then 


P+(x0,k) 


<t>+{xo,k)  _  1  -  p+(k)e  2,kxo 
ik<f>+(x0,k)  l -\- p+(k)e~2'kx° 

gikx  0 

<f>+(xo,  k) 


(2.30) 
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Similarly,  for  any  Xi  >  1, 


P-{x  1,*) 


1  -  fi-(k)e2ikxi  _ 

1  +  fi-(k)e2lkxi  4>-(xi,k) 


(2.31) 


Chapter  3 

Mathematical  Preliminaries 


In  this  chapter,  we  summarize  several  well-known  mathematical  facts  to  be  used 
in  the  rest  of  this  thesis.  These  facts  are  given  without  proofs,  since  Lemmas  3.1- 
3.8  are  found  in  standard  textbooks  (see,  for  example,  [10],  [11],  )  and  Lemmas 
3. 9-3. 6  are  easy  to  verify  directly. 

3.1  Basic  Lemmas 

The  following  basic  facts  are  tailored  and  stated  in  such  a  way  that  they  will  be 
directly  used  in  the  existence  and  convergence  proofs,  see  Chapter  4  and  Section 
5.2. 

Lemma  3.1  Suppose  that  A  is  a  linear  mapping  C[0, 1]  — ►  C[0, 1]  and  that 
ll-AH  <  p,  with  p  a  real  number  such  that  p  <  1.  Then  for  any  g  €  L2[0, 1], 
the  equation 

<f>  =  A<j>  +  g  (3.1) 

has  a  unique  solution,  which  is  the  sum  of  the  series  (known  as  Neumann’s  series) 

00 

*  =  Y,A>S-  (3.2) 

j=0 

Furthermore, 

n  n+l 

U-Y,A’9\\<(—  M-  (3.3) 

Lemma  3.2  Suppose  that  f  €  c^([0,  Z)])  (i.e.,  f  has  m  continuous  derivatives 
and  f(x)  =  0  for  all  x  £  (0,-0)^,  and  that  /h")  is  absolutely  continuous.  Suppose 
further  that  g  €  cm+1(R),  is  absolutely  continuous  and  there  exist  real 
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number  a  >  0  such  that  g'(x )  >  a  for  all  x  €  R.  Then  there  exists  a  real  c  >  0 
such  that 

(3.4) 

for  all  complex  k  such  that  Im(k )  >  0. 

Lemma  3.3  Suppose  that  f  €  cl(R )  with  l  a  nonnegative  integer.  Suppose  fur¬ 
ther  that  fW( 0)  =  0  for  0  <  j  <  l,  /O)  is  absolutely  continuous  Then  there 
exists  a  positive  number  c  such  that 

/o/«)e“<^  =  -g(^)  /'-(*)+ (gE)  (/(*)  +  •(«.*))  (3-5) 

6  :  R  x  C+  — >  C  an  absolutely  continuous  function  of  x  €  [0, 1]  sticft  that 

\b(x,k)\<c.  (3.6) 

/or  all  x  €  [0, 1],  k  €  C+ .  Furthermore,  if  f(x)  =  0  for  all  x  >  D  with  D  a 
positive  number,  then 

|6(a;,  k)\  <  c.  (3.7) 

for  all  (x,  k)  6  R  x  C+ 

Lemma  3.4  Suppose  that  a  :  [0, 1]  R  and  b  :  [0, 1]  — *  C  are  two  continuous 
functions,  and  that  a(x)  >  0,  for  all  x  €  [0, 1].  Then  for  any  two  solutions  u  and 
v  of  the  second  order  ODE 

(a(x)<t>'(x))'  +  b(x)<f>{x)  =  0,  (3.8) 

there  exists  a  constant  c  such  that 

a(x)(u(x)v'(x)  —  v(x)u'(x))  =  c  (3.9) 

for  all  x  €  [0, 1] .  Furthermore,  c  ^  0  if  and  only  if  u  and  v  are  linearly  inde¬ 
pendent.  (The  expression  W(u,v )  =  u( x)v'(x)  —  u(x)u'(x)  is  referred  to  as  the 
Wronskian  of  the  pair  u,v). 

Lemma  3.5  (Gronwall’s  inequality)  Suppose  that  u,v,w  :  [0,a]  — >  R  are  three 
continuous  and  nonnegative  functions,  satisfying  the  inequality 

w(x)  <  u(x)  +  f  v(t)w(t)dt  (3.10) 

Jo 

for  all  x  €  [0,  a] .  Then 

w(x)<u(x)+  I  u(t)v(t)ef'  v^drdt  (3.11) 

JO 


for  all  x  €  [0,a]. 
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The  following  lemma  is  a  special  case  of  the  general  theorem  about  continuous 
dependence  on  initial  conditions  and  parameters  of  solutions  of  ODEs  (see,  for 
example,  [11]). 

Lemma  3.6  Suppose  that  a  :  C  — ►  C  is  an  entire  function  and  that  A  :  RxC  — ► 
CnXn  is  an  n  x  n-matrix  whose  entries  aiffx,  z),i,j  =  1 ,  ...,n  are  continuous 
functions  of  x  and  entire  functions  of  z  for  all  x  €  R.  Then  for  any  z  €  C,  the 
differential  equation 

y'(x,z)  =  A(x,z)  ■  y(x,z)  (3.12) 

subject  to  the  initial  condition 

y{  0)  =  c(z)  (3.13) 

has  an  unique  solution  y(x,z )  for  all  x  €  R.  Moreover,  y(x,z)  is  an  entire 
function  of  z. 


3.2  Schrodinger  Equation  and  Riccati  Equa¬ 
tion 

Lemmas  3.7-  3.13  describe  the  basic  facts  about  the  Helmholtz  equation  and 
its  connections  with  the  Schrodinger  Equation  and  the  Riccati  Equation,  in  the 
context  of  scattering  problems. 

Lemma  3.7  Suppose  that  Gk  :  [0, 1]  x  [0, 1]  — *•  C  is  the  Green’s  function  of  the 


boundary  value  problem 

i>"(x,  k)  4-  k2xf(x,  k)  =  0,  (3-14) 

xl>'(Q,k) ikil>(0.,k)  —  0,  (3.15) 

—  iki(>(l,k)  =  0.  (3.16) 

for  any  complex  k  ^  0.  Then  the  boundary  value  problem 

tf"{x,  k )  +  (k2  +  ?7(x))^(x,  k )  =  /(x,  k )  (3.17) 

^'(0,  k)  +  ikxl>(0,  k)  =  0,  (3.18) 

rp'(l,k)  —  ikth(l,k)  =  0.  (3.19) 

is  equivalent  to  a  second  kind  integral  equation 

tj>(x,k)  =  -f  Gk{x,t)r)(t)ip(t,k)dt  +  g{x,k)  (3.20) 

Jo 


with  f,g  :  [0, 1]  x  C  — >  C  and  g  defined  by  the  formula 

g{x,k)=f  Gk{x,t)f(t,k)dt. 
Jo 


(3.21) 
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Lemma  3.8  For  any  complex  k  0,  the  Helmholtz  equation 

ip"(x,  k )  4-  k2ip(x ,  k)  =  0  (3.22) 

with  the  outgoing  radiation  conditions  (2.8)  (2.9)  has  the  Green’s  function 


1  (  x  <  t 

^ x  i !: 


(3.23) 


Lemma  3.9  Suppose  that  q  :  R  — >  R  is  a  c2— function  such  that  q  >  —  1  for  all 
x  E  R.  Suppose  further  that  the  functions  n,x,S,ri,g  :  R  — *  R  are  defined  by  the 
formulae 

(3.24) 


v(t)  = 


n(x)  =  y/l  +q(x), 
t(x)  =  /  n(r)c/r, 

S(t)  =  (1  +  q(x(t))r< 
S"{t)  n'(x) 


(3.25) 

(3.26) 


S(t)  2(n(x))2 

=  ^(i  +  q)~2  (q"(x)  -  \/i  +  q(xW{x)  -  ^(i  +  9)~V2(*))  (3-27) 


gif)  =  =  f(x)  •  (1  +  $(*))* 

Finally,  suppose  that  the  function  <f>  :  R  x  C  — »  C  satisfies  the  equation 
cf)"(x,  k )  +  fc2(l  +  ?(x))  •  <f>(x,  k )  =  /(x), 
and  the  function  ip  :  R  x  C  — >  C  is  defined  by  the  formula 

ip{t,  k )  =  <p{x(t),  k)/S(t )  =  <p(x,  k)  •  (1  +  g(x))< 

Then  the  function  ip  satisfies  the  Schrodinger  equation 

ip"(t,  k )  +  ( k 2  +  rj(t))  •  ip(t,  k )  =  g(t). 


(3.28) 


(3.29) 


(3.30) 


(3.31) 


at  all  t  €  R. 

Remark  3.10  Lemma  3.9  provides  a  connection  between  the  solutions  of  the 
Helmholtz  equation  (3.29)  and  those  of  the  appropriately  chosen  Schrodinger 
equation  (3.31).  This  connection  will  be  used  in  the  following  chapter  as  an 
analytical  tool.  However,  it  is  not  useful  in  numerical  computations  since  the 
connection  between  rj  and  q  (see  (3.27))  is  generally  ill-conditioned. 
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Corollary  3.11  Suppose  that  under  the  conditions  of  the  preceding  lemma  that 
q(x)  =  0  for  all  x  £  (0, 1).  Suppose  further  that  the  functions  V>+,  xb.  :  RxC  — »  C 
are  defined  by  the  formulae 


xl>+(t,k)  =  <j>+(x(t),k)/S(t), 

(3.32) 

xf-{t,k)  =  <j>-(x(t),k)/S(t). 

(3.33) 

Then  xp+,  xp-  satisfy  the  ODEs 

xp'l(t,k)  +  (k2  +7}{t))  -x/>+(t,k)  =  0, 

(3.34) 

xp'L{t,  k)  +  (k2  +  r](t))  ■  xf_(t,  k)  =  0 

(3.35) 

subject  to  the  boundary  conditions 

xl>+(t,k)  =  tik)-eik{t-Tl) 

(3.36) 

for  all  t  >T\,  and 

xP-{t,k)  =  e~ikt 

(3.37) 

for  allt  <  0  with  T\  >  0,  ffik)  0  defined  by  the  formulae 

T\  =  t(l)  =  f  n(r)dr, 

JO 

(3.38) 

{(i)  =  S(T,)eik. 

(3.39) 

Furthermore, 

_  /_  ;r\  _  r/r>  V>V(<,  k)  n'{x) 

p+(x,k)  -  n{*)ik^k)  2 ikn(Xy 

(3.40) 

p-(x,k)  =  n(x)~  , 

~ikxf_(t,  k)  2 ikn(x) 

(3.41) 

Observation  3.12  Suppose  that  q(x)  =  0  for  all  x  (0,1).  Then  according  to 
Lemma  3.9  and  Corollary  3.11, 


t  =  x, 

S(t)  =  1, 


(3.42) 

(3.43) 


and  consequently 


4>+(x,k)  =  x l>+{t,k) 


for  all  x  <  0.  Now,  suppose  the  function  x fi+  is  defined  by  formulae 
(2.1 7).  Defining  the  scattered  field  Vw+  :  R  x  C  — ►  C  by  the  formula 


(3.44) 

(3.32), 


xp+(t,  k)  =  e'kt  +  Vw+(f>  k). 


(3.45) 
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we  immediately  see  that 

i’scat+it,  k)  =  fi+(k)  ■  e~lkt  (3.46) 

for  all  x  <  0  due  to  (3-44),  (3-17),  (3.45).  Finally,  combining  (3-45)  with 
(3.34),  we  observe  that  ^5Cat+  satisfies  the  Schrodinger  equation 

—  k2a( 

Cal+(*>  k)  +  (*2  +  y(t)W,cat+(t,  k)  = - 2LL2J -  (3.47) 

subject  to  outgoing  radiation  conditions  (2.8),  (2.9)  (the  latter  due  to  (3-46), 
(3.36)). 

Lemma  3.13  Suppose  that  under  the  conditions  of  the  preceding  lemma, 


4>+(x0,ko)  ±  0,  (3.48) 

<j>-(x0,k0)  ±  0  (3.49) 

at  some  point  (x0,k0)  €  R  x  C .  Then  there  exists  a  neighborhood  D  of  (x0,ko) 
such  that  the  impedance  functions  p+,p~  satisfy  the  Riccati  equations 

p'+(x,k )  =  -ik(p2+  (x,  k)  -  (1  +  g(x)),  (3.50) 

p'_(x,k)  =  ik(p2_ (x,  k)  —  (1  +  q(x))  (3.51) 

for  all  (x,  k )  €  D. 

Observation  3.14  Combining  formulae  (2.23),  (2.24),  we  easily  observe  that 

p+(x,k )  =  1,  for  all  x  >  1,  (3.52) 

p_(x,fc)  =  1,  for  all  x  <  0,  (3.53) 


for  all  complex  k  /  0. 


Chapter  4 

Impedance  Functions  and  Their 
Properties 


In  this  chapter,  we  investigate  analytical  properties  of  the  impedance  functions 
p+ , p_ .  Our  principal  purpose  here  is  to  formulate  exactly  and  prove  the  following 
three  facts. 

(1)  For  any  x  €  R,  the  impedance  functions  p+(x,  fc),p_(x,  k)  are  analytic  func¬ 
tions  of  k  in  the  upper  half  plane  C+.  Furthermore, 

p+(x,k)  =  +  +  (4.1) 

p.(x,k)  =  l/l +?(*)+ 4  (l+^(l))-^  +  CI(t-2).  <4-2) 

for  all  x  €  R,  k  €  C+  (see  Theorem  4.14  below). 

(2)  For  large  real  k,  the  difference  between  p+  and  p_  is  extremely  small  (it 
decays  like  k~m ,  where  m  is  the  smoothness  of  the  scatterer,  see  Theorem  4.18 
below).  The  expressions  (4.1),  (4.2)  are  the  first  two  terms  in  WKB  expansions 
of  the  functions  p+ ,  p_ ,  respectively. 

(3)  For  any  a  >  0,  and  all  x  G  R,  we  have  the  so-called  trace  formula 

?'(*)  =  §(!  +  ?(*))  /  (p+(*,*)  -p-(x,k))dk  +  0{a-{m-1}),  (4.3) 

IT  J-a 

with  m  the  smoothness  of  the  scatterer  (see  Theorem  4.19  below). 

As  often  happens,  the  statements  (l)-(3)  above  have  extremely  simple  for¬ 
mulations,  and  a  transparent  physical  interpretation.  However,  their  proofs  are 
technical  and  do  not  follow  any  simple  physical  intuition. 
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4.1  Boundedness 

The  following  five  lemmas  establish  the  basic  properties  of  the  impedance  func¬ 
tions  p+,p~  introduced  in  Chapter  1.  Lemma  4.1  is  a  technical  one,  describing 
the  behavior  of  <f>+,  4>-  in  the  vicinity  of  k  =  0  in  the  complex  plane.  Lemma  4.2 
describes  the  properties  of  the  impedance  functions  p+,p~  near  A:  =  0,  Lemma 
4.4  demonstrates  the  well-definedness  of  the  impedance  functions  for  real  k,  and 
Lemmas  4.5  and  4.6  provide  upper  and  lower  bounds  for  the  impedance  func¬ 
tions. 

Lemma  4.1  Suppose  that  q  £  c([0,l])  and  A  >  0  is  a  real  number.  Then  there 
exist  three  positive  numbers  6,  a  and  0  such  that 


1. 

\<f>+{x,k)  -  1|  <  a|fc|, 

(4.4) 

2. 

\<f)-(x,k)  -  1|  <  a|fc|, 

(4.5) 

3. 

\<f>'+(x,k)-ik\<0\k\2. 

(4.6) 

4. 

\<f>'_(x,k)  +  ik\<m2, 

(4.7) 

5. 

(j>+(x,k)  ±  0, 

(4.8) 

6. 

d)-(x,  k )  ^  0, 

(4.9) 

for  all  real  x  £  [—A,  A]  and  complex  k  such  that  j&|  <  S. 

Proof.  Since  the  proofs  of  this  lemma  for  <f>+,  <j>'+  and  for  cf>_,  <j>'_  are  identical, 


we  only  prove  it  in  the  case  of  </>_,  <p'_.  Defining  two  auxiliary  functions  : 

R  x  C  — *  C  by  the  formulae 

<t>i{x,k)  =  (j>~(x,  k)  —  1,  (4.10) 

xp(x.k)  =  <f>’_(x.  k)  +  ik,  (4.11) 

and  combining  (4.10),  (4.11)  with  equation  (2.1)  and  the  initial  condition  (2.27), 
we  observe  thac  the  functions  <f> satisfy  the  linear  first  order  ODEs 

^[(x^k)  =  tl>(x,k)  +  ik ,  (4.12) 

4>'(x,  k)  =  -k2(l+q(x))(l  + (j>i(x,k))  (4.13) 

subject  to  the  initial  conditions 

&(0,*)  =  0,  (4.14) 

^(0,  k)  =  0.  (4.15) 


We  start  with  showing  that  there  exist  continuous  functions  M,  N  :  RA  x  RA  — » 
R+  such  that,  for  any  s  £  R+,  M(s,t),  iV(,s,f)  are  monotonicaily  increasing 
functions  of  t  for  all  t  £  R+  and 

<  M(A,\k\)\k\, 

hK*>*)l  £  jv(a,  |fc|)|fc|2. 


(4.16) 

(4.17) 
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First,  we  prove  the  estimate  (4.17).  Integrating  (4.12)  from  0  to  x,  we  have 

<f>i(x,k)=  f  (ik  +  ip(t,k))dt,  (4-18) 

Jo 

and  substituting  (4.18)  into  (4.13)  and  integrating  the  result  of  the  substitution, 
obtain 

V>(x,  k)  =  - k 2  jf  (1  +  q(t))  (l  +  J  (ik  +  VK'D  k))dr^  dt.  (4.19) 

Denoting  \x/>(x,k)\  by  a(x,k )  and  observing  that  1  +  q(x)  <  n\  (see  (2.3)  in 
Section  2.1),  we  obtain 

a(x,k)  <  {k^nl  ^\x\  +  ^x2\k\  +  J  J  a(r,  k)drdt ^ 

<  \k\2nl(^\x\  +  ^x2\k\j  +  \k\2nl  J  (x  -  t)a(t:k)dt  (4.20) 

for  any  x  €  R-  Gronwall’s  inequality  (see  Lemma  3.5)  implies  that  for  any 
x  €  [0,A], 

a(x,k)  <  \k\2n\(\x\  +  \>x2\k\+  \t\+ ^t2\k\(x  -  t)e^(x~t)2 dt^ 

<  N(A,  |fc|)|fc|2.  (4.21) 

It  is  easy  to  see  that  (4.21)  is  also  valid  for  any  x  €  [ — j4,  0],  and  we  obtain  the 
estimate  (4.17)  with  N(A ,  k)  defined  by  the  formula 

N(A,  |*|)  =  sup  nl(\x\  +  ~x2\k\+  f  |*|  +  ^2|&|(a:  -  tje&'^dt)  .  (4.22) 

-A<x<A  \  Z  J 0  2  / 


We  now  turn  our  attention  to  the  estimate  (4.16).  Substituting  (4.17)  into 
(4.18),  we  obtain 


l*(*,*)l  <  M  (1*1  +  N(A,  |*|)|*f) 

<  M(A,\k\)\k\,  (4.23) 

with 

M(A,  |*|)  =  A(1  +  |*|JV(A,  |*|)),  (4.24) 

for  all  real  x  6  [—A,  A]  and  complex  k,  which  proves  (4.16). 

Now,  the  estimates  (4.5)  and  (4.7)  easily  follow  from  (4.16)  and  (4.17). 
Indeed,  since  M(A,t)  is  a  continuous,  monotonically  increasing  function  of  f, 
there  exists  a  real  6  such  that 


M(A,S)-6<  1. 


(4.25) 
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Denoting  M(A,8)  by  a ,  N(A,8)  by  0  and  observing  that  M(A,  |&|),  N(A,  |fc|) 
are  monotonically  increasing  functions  of  \k\,  we  have 

|^i(*,Jb)|  <  M{A,\k\)\k\<M(A,8)\k\  =  a\k\,  (4.26) 

W*,*)l  <  ml*l)|i|<AW)l*l  =  /»l*|J,  (4.27) 

from  which  (4.5),  (4.7)  follow  immediately. 

Finally,  (4.9)  is  a  direct  consequence  of  (4.26)  and  (4.25).  □ 

Lemma  4.2  Suppose  that  q  €  Cq([0,  1])  and  A  >  0  is  a  real  number.  Then  there 
exists  8  >  0  such  that  the  impedance  functions  p+ .  p_  are  continuous  functions 
of  (x,  k)  for  all  real  (x,  k)  £  D  with 

D  =  {(x,  Jfe)|x  €  [-A,  A\,keC,k  ±  0,  |Jb|  <  <5}  (4.28) 

Furthermore, 

limp+(x,A:)  =  1,  (4.29) 

limp_(x,&)  =  1.  (4.30) 

Proof.  Due  to  Lemma  4.1,  there  exists  a  positive  number  8  such  that 

<£+(x,  k)  7^  0,  <f>-(x,  k)  ^  0  for  all  real  (x,  k)  €  D.  Therefore,  the  functions  p+,p~ 
are  well-defined  in  D,  and  their  continuity  follows  from  the  continuity  of  <f>+, 
<f>+,  <t>-,  <t>'~  and  the  formulae  (2.15),  (2.16).  Finally,  (4.29),  (4.30)  are  direct 
consequences  of  Formulae  (4.4)-  (4.7).  □ 

Remark  4.3  While  the  impedance  functions  p+,p~  are  continuous  in  the  vicinity 
of  k  =  U  in  the  complex  plane,  formulae  (2.15),  (2.16)  fail  to  define  p+,p~  at 
k  =  0.  We  now  can  define  p+(x,0)  =  p_(x,0)  =  1  for  all  x  €  R  due  to  Lemma 
4.2. 


Lemma  4.4  For  any  real  k  /  0  and  all  x  €  R 

d>+(x,k)  ^  o, 

(4.31) 

<t>'+(x,k)  ±  0, 

(4.32) 

<t>-{x,k)  ^  0, 

(4.33) 

<f>'_(x,k)  ^  0. 

(4.34) 

Proof.  Again,  since  the  proofs  of  this  lemma  for  <f>+,<j>'+  and  for  d>_,  <f>'_  are 
identical,  we  only  prove  (4.33)  and  (4.34).  Denoting  the  real  part  of  <5-  by  u 
and  the  imaginary  part  by  v,  so  that 


<f>-(x,k)  =  u(x,  k)  +  zt>(x,  k), 
<t>'_(x,k)  =  it'(x,  k)  +  iv'(x,  k), 


(4.35) 

(4.36) 
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we  observe  that  each  of  the  functions  u,v  satisfies  equation  (2.1)  (since  the 
coefficients  of  the  equation  are  real).  Combining  the  initial  condition  (2.27) 
with  (4.35),  we  immediately  see  that 

u(x,k)  =  cos  (kx),  (4-37) 

v(x,  k)  =  sin(kx)  (4.38) 

for  all  x  <  0  and  k  ^  0.  Therefore,  the  Wronskian  of  the  pair  u,  v  is 

W(u,v)  =  k,  (4.39) 

for  any  x  €  R  (see  Lemma  3.4),  and  u(x,k),  v(x,  k)  can  not  be  both  zero,  nor 
can  u\x ,  k),  v'(x,  k),  for  any  x  €  R  and  k  ^  0.  Now,  formulae  (4.33)  and  (4.34) 
immediately  follow  from  (4.35)  and  (4.36)  □ 

We  have  shown  that  the  impedance  functions  p+,p_  are  well-defined  for  all 
real  k  (see  Lemmas  4.2,  4.4  and  Remark  (4.3)).  Now,  we  turn  our  attention  to 
the  well-definedness  of  the  impedance  functions  on  the  upper  half  of  the  lb-plane. 
First  we  provide  the  lower  bounds  for  p+,p_. 

Lemma  4.5  For  all  x  €  R  and  any  k  such  that  lm(k)  >  0, 

Re(p+(x,k))  >  n0sin(arg(£)),  (4.40) 

i?e(p_(x,  k))  >  n0sin(arg(fc))  (4.41) 

with  0  <  n0  <  1  the  minimum  of  n(x)  ( see  (2.3)  in  Section  2.1),  and  arg(fc) 
the  argument  of  the  complex  wave  number  k. 

Proof.  Since  the  proof  of  (4.40)  and  that  of  (4.41)  are  identical,  we  only 
provide  the  latter.  Observing  that 

Re(p_(x,k))  =  p_(x,fc)  =  1  >  rig sin(arg(fc)),  (4.42) 

for  any  Im(k )  >  0  and  all  x  <  0  (see  (3.53)  in  Chapter  3),  we  will  prove  (4.41) 
by  showing  that 

~  (i?e(p_(x,  k)))  >  0  (4.43) 

for  any  x  >  0  such  that 

0  <  i?e(p_(x,  k))  <  n0  sin(arg(fc))  (4.44) 

(obviously,  0  <  arg(k)  <  v  for  any  k  such  that  Im(k)  >  0). 

We  will  denote  by  a,  6,  it,  v  the  real  and  imaginary  parts  of  k  and  p_  respec¬ 
tively,  so  that 

k  =  a  +  ib, 

p_(x,fc)  =  u(x,  k)  +  iv(x,k), 


(4.45) 

(4.46) 
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with  b  >  0.  Now,  we  can  rewrite  the  Riccati  equation  (3.51)  for  p_  in  the  form 

u'  —  b(v2  —  u2  +  n2)  —  2auv,  (4.47) 

v  —  —a(v2  —  u2  +  n 2)  —  2  buv.  (4.48) 

We  observe  that  k )  is  a  function  of  u,  v  given  by  the  formula 


Q 

—u(x,  k)  —  /(u,  u)  =  b(v2  —  u2  +  n 2)  —  2auv.  (4.49) 

ox 

Denoting  the  interval  [0,  n0  sin(arg(fc))]  by  /,  and  defining  the  region  D  C  R  x  R 
via  the  formula 

D  =  {(u,t>)|u  €  I,v  e  i?},  (4.50) 

we  observe  that 

min  f(u,v )  =  b(n 2  —  rig)  >  0  (4.51) 

(u,ti  )€D 

which  proves  (4.43)  given  (4.44).  Now,  (4.41)  follows  immediately  form  (4.42), 
(4.43)  and  (4.44).  □ 

As  a  direct  consequence  of  Lemma  4.5,  the  following  lemma  establishes  the 
upper  bounds  of  the  impedance  functions  in  the  upper  half  plane. 

Lemma  4.6  For  any  k  such  that  Im(k )  >  0  and  all  x  €  R, 


n  i 


[(?+(*.  *))l  <  sin(arg (i))’ 

l(p-(*,fc))l  <  .  ,n',.Tv 

sin(arg(fc)) 

with  n\  >  0  the  maximum  ofn(x)  (see  (2.3)  in  Section  2.1). 

Proof.  Again,  we  only  give  the  proof  of  (4.53)  since  the  proof  for 
identical.  According  to  Lemma  4.5,  the  function 


(4.52) 

(4.53) 


(4.52)  is 


r(x,k)  =  l/p„(x,k)  (4.54) 

is  well-defined  for  any  Im(k)  >  0.  Combining  (4.54)  with  the  equation  (3.51) 
and  the  boundary  condition  (3.53)  for  p_,  we  observe  that  r(x,  k)  obeys  the 
Riccati  equation 

r'(x,  k )  =  ikn2{x)  ^r2(z,  fc)  -  j  ■,  (4.55) 

subject  to  the  initial  condition  r(0,  k)  =  1.  Reproducing  the  proof  of  Lemma  4.5 
almost  verbatim,  we  obtain  a  lower  bound  for  the  real  part  of  r 

r(x.fc))  > 

Tli 


(4.56) 
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|(P.(*,^)l<fie(r(*,*))  —  sin(arg(fc)) 

(4.57) 

is  readily  obtained  by  combining  (4.54)  with  (4.56).  □ 

Corollary  4.7  For  all  x  €  R  and  k  such  that  Im(k)  >  0, 

<t>+(x,k)  ±  0, 

(4.58) 

4>+(x,k)  ^  0, 

(4.59) 

<t>-{x,k)  ±  0, 

(4.60) 

<f>'_(x,k)  ±  0. 

(4.61) 

Proof.  We  prove  this  corollary  by  contradiction.  First,  we  observe  that 

<p+(x,  k)  =  0 

(4.62) 

implies 

<f>+(x,  k)  =  0 

(4.63) 

and  vice  versa,  since  both  <f>+(x,k)  and  <f>'.(x,  k)  are  continuous  functions  of  x, 

and  their  ratio 

‘*-p+(*M=  Mxk) 

(4.64) 

is  bounded  from  both  above  and  below  due  to  Lemmas  4.5, 

4.6. 

Suppose  now  that  for  some  x0  €  R,  Im(k0)  >  0, 

<t>+(xo i  k0)  =  4>'+{x0,  ho)  =  0. 

(4.65) 

Then  the  pair  of  functions 

4>{x)  =  (j>+(x,k0), 

(4.66) 

ip{x)  =  4>'+{x,k0) 

(4.67) 

satisfies  the  system  of  ODEs 

<j>'(x)  =  0(x), 

(4.68) 

^'(x)  =  -^oC1  +  ?(*))<£(*), 

(4.69) 

subject  to  the  initial  conditions 

<t>{x0)  =  xj}(x0)  =  0. 

(4.70) 

However,  the  initial  value  problem  (4.68),  (4.69),  (4.70)  has  a  unique  solution 


(j>(x)  =  x/;(x)  =  0  (4.71) 

for  all  x  <  x0,  which  contradicts  the  condition  7 2.26),  proving  (4.58),  (4.59). 
The  proof  of  (4.60),  (4.61)  is  identical.  O 
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Observation  4.8  Due  to  Lemma  3.13,  it  is  easy  to  see  that 

P+(x »&)  =  P+  (*,-£),  (4.72) 

p-(x,k)  =  p_(x,-T),  (4.73) 

for  all  x  G  R  and  k  G  C+ .  For  real  k,  equalities  (4-72),  (4-73)  assume  the  form 

P+(x,k)  =  P+(xi~k),  (4-74) 

p~(x,  k)  =  p-(x,-k).  (4.75) 

Indeed,  combining  the  complex  conjugate  of  { 3.50 1  with  that  of  (3.52),  we  obtain 
the  ODE 

(p+(x,fc))  -  -i(-k)(p+(x,kf  -  (1  +  ?(*))  (4.76) 

subject  to  initial  condition 

P+(0,k)  =  1.  (4.77) 

Now,  replacing  k  by  — in  (3.50)  znd  (3.52),  we  have 

p'¥(x, -k)  =  -*'(-t)(p+ («,  - k )2  -  (1  +  q(x)),  (4.78) 

and 

p+(0,-k)  =  l.  (4.79) 

We  notice  that  p+(x,k),  p+(x,  —k)  satisfy  identical  differential  equations  (4-76), 
(4-78)  with  identical  boundary  conditions  (4- 77),  (4-79),  from  which  (4-72) 

follows.  A  similar  calculation  proves  (4-73). 

4.2  Smoothness  and  Asymptotics 

The  following  two  technical  lemmas  describe  the  asymptotic  behavior  of  the  func¬ 
tions  ip+,xh-  (see  Corollary  3.11  in  Chapter  3),  4>+  and  <f>-  at  large  frequencies. 
They  will  be  used  in  proofs  of  Theorems  4.14,  4.18,  describing  the  high-frequency 
asymptotics  of  the  impedance  functions  p+ ,  p-.  Theorems  4.14,  4.18  are  in  turn 
used  in  the  following  chapter  to  derive  the  trace  formulae  (4.194),  (4.198),  which 
are  the  principal  analytical  tool  of  this  thesis. 

Lemma  4.9  Suppose  that  for  any  a  >  0,  the  region  K(a)  C  C  is  defined  by  the 
formulae 

K(a )  =  €  C ,  Im(k )  >  0,  |&|  >  a}.  (4.80) 

Suppose  further  that  q  G  Cq([0,  1]),  q(x)  >  —1  for  all  x  G  R,  and  the  second 
derivative  of  q  is  absolutely  continuous.  Then  there  exist  real  numbers  A  > 
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0,  c  >  0  such  that 

*+(«,*)  =  «fc)ei‘"-T'>(l  +  ^_(\(T)<iT  +  e+(<,J:)),  (4.81) 

rl>'+(t,k)  =  ik((k)e‘^-T'l(l  +  ~j'rilT)dr+6+(t,k)y  (4.82) 

4i-(t,k)  =  e-i“  (l  +  l(^)dT  +  £-(<.*))  ,  (4.83) 

<f'Jt,k)  =  -ike-iU(l  +  ±J\(T)dT  +  64t,k)y  (4.84) 

with  £(k)  :  C  — *  C,  T\  >  0  defined  by  (3.39),  (3.38)  (see  Corollary  3.11  in 
Chapter  3),  and  e+,e_,6+,6_  :  R  x  K(A)  — >  C  continuous  functions  such  that 

|e+(t,fc)|  <  c-k~2,  (4.85) 

\6+(t,k)\  <  c  •  k~2,  (4.86) 

\e.(t,k)\  <  c  •  k~2,  (4.87) 

\S-(t,k)\  <  c-k~2.  (4.88) 

for  all  (t,  Jc)gRx  K(A). 

Proof.  Since  the  proofs  of  this  lemma  for  ip'+  and  for  are  identical, 

we  only  prove  it  in  the  case  of  t/*-,  V’L-  Introducing  two  auxiliary  functions  m,  n  : 
R  x  C  — >  C  by  the  formulae 

m(t,k)  =  elktil>-(t,k),  (4.89) 

n(t,k)  =  -—eikti(>'_{t,k)  (4.90) 

ik 

and  combining  (3.35),  (3.37)  with  (4.89),  (4.90),  we  observe  that  m  satisfies 
the  equation 

m"(t,  k )  —  2 ikm'(t,  k)  =  —r}(t)m(t,  k)  (4-91) 

(77  £  co([0,Ti])  is  absolutely  continuous,  see  (3.27)  for  the  definition  of  77)  subject 
to  the  initial  conditions 

m(0,  k)  —  1,  (4.92) 

m'{0,k)  =  0.  (4.93) 

Multiplying  (4.91)  by  e~2lkt  and  integrating  the  result  from  0  to  t,  we  have 

m’(t,k)  =  —  f  Tj(T)e2,k^~r^m(T,k)dT. 

Jo 


(4.94) 
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Integrating  (4.94)  from  0  to  t,  we  obtain  the  second  kind  Volterra  integral  equa¬ 
tion  for  m 

m  =  Fk(m)  +  1  (4.95) 

with  the  mapping  F*  :  c(R )  — »  c(R )  defined  by 

A(/)(i)  =  ^j‘  xM(l  -  e2i«-r'mT)dT.  (4.96) 

Combining  (4.94)  with  (4.89),  (4.90),  we  observe  that 

n(t,k)  =  m(t,k)  —  —r  f  q{T)e2tk^~T^m{T,k)dT.  (4.97) 

2z  k  J o 

Since  77  €  co([0,  Ti]),  the  function  t?(t)(1  —  e2,k^~T^)  is  bounded  for  all  real  t,r 
and  k  €  F(0).  Therefore,  there  exists  a  real  number  c\  >  0  such  that 

II All  <  j|i ,  (4-98) 

and  hence  there  exists  a  real  number  A  >  0  such  that 

HAH  <  1  (4.99) 

for  all  k  €  K{A).  Now,  according  to  Lemma  3.1,  for  all  ( t,k )  €  R  x  K(A ), 
the  unique  solution  of  (4.95)  can  be  approximated  by  the  Neumann’s  series 
truncated  at  the  second  term 

m(t,k)  =  l  +  ijf\(T)(l-e2i“'-'>)<iT  +  Q(i,t) 

=  1  +  j(  X(A*  +  0(t,  k)  +  o(t,k)  (4.100) 

with  a,  0  :  Rx  K(A)  — *  C  such  that 

IMI  <  ^  (4.101) 

(see  Lemma  3.1),  and 

0(t,k)  =  -- \lHn(TVik{'-')dT-  (4.102) 

Since  q"  is  absolutely  continuous  and  q(x)  =  0  for  all  x  <  0,  we  observe  that  77 
is  absolutely  continuous  and  77(2:)  =  0  for  all  x  <  0  (see  (3.27)  in  Lemma  3.9). 
According  to  Lemma  3.3,  there  exists  C2  >  0  such  that 
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for  all  x  £  [0,1],  k  £  C+.  Now,  combining  (4.100)  with  (4.101)  and  (4.103), 
we  observe  that  there  exists  c$  >  0  such  that 


m(U)"(1  +  2 IkSl^7) 


-  |*|* 


for  all  ( t ,  k)  £  R  x  /f(v4).  Similarly,  there  exists  c4  >  0  such  that 


<  -ii- 
-  |*|’ 


(4.104) 


(4.105) 


due  to  (4.97),  (4.104). 

Now,  (4.83),  (4.87)  follow  immediately  from  (4.104),  (4.92),  and  (4.84), 
(4.88)  are  a  direct  consequence  of  (4.105),  (4.93).  □ 


Lemma  4.10  Suppose  that  q  £  Cq([0,  1]),  7  >  2,  is  absolutely  continuous 
and  q{x)  >  —1  for  all  x  £  R.  Then  for  any  integer  1  <  l  <  7,  the  l-th  interate 
mi  :  R  x  C+  —*  C  defined  by  the  formulae 


m0(t,k)  =  0,  (4.106) 

m/(t,fc)  =  1  +  Fk(mi-i)(t,k)  (4.107) 

=  1  +  2 lie  Jo  ~  k)dr  (4.108) 

(see  (4-95),  (4-96))  assumes  the  form 

mi(t,  k)  =  1  +  £  aj(t)  +  °7(t,  k)  (4.109) 


with  aj  :  R  —+  R,  j  =  1, . . .  ,7  —  1,  a7  :  R  x  C+  — *  C  such  that 

dx^ 


(4.110) 


are  bounded  and  absolutely  continuous  for  all  x  £  R,  j  —  1, ...  ,7  —  1,  and 


<h(t,k)  (4.111) 

is  bounded  and  absolutely  continuous  function  oft  for  all  (t,k)  £  R  x  C+. 

Proof.  We  prove  this  lemma  by  induction.  For  1=1,  formulae  (4.106), 
(4.108)  yield 

mi(t,fc)  =  l  (4.112) 

for  all  ( t ,  k)  £  Rx  C+ ,  which  is  already  in  the  form  (4.109)  satisfying  conditions 
(4.110),  (4.111). 
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For  /  >  1,  assuming  that  is  in  the  form  (4.109)  satisfying  conditions 

(4.110),  (4.111),  we  obtain  m;+i  using  (4.108): 

m/+i(t,  k)  =  1  +  2k  i)  7?(T)(1  “  e2'k(t~T])mi{T,k)dT 


—  1  +  I\(t ,  k )  -+■  /2(t,  k )  +  Iz(t,  fc)  +  I${t,  k) 

(4.113) 

with  Ij  :  R  x  C+ 

— : >  C,  1  <  j  <  4  defined  by  the  formulae 

h(t,k )  = 

2k  /o’ ,(T)<iT  +  g  (2k)1  /„’ 

(4.114) 

h(t,k)  = 

-Li  *»-**->)*. 

(4.115) 

h(t,k)  = 

(^)  jf  ??(r)as-i(T)e2,fc(t  T)dr  =  -  J]  Js(f,  *), (4.116) 

U (t,  k )  = 

(4.117) 

Clearly,  we  only  need  to  show  that  Ij ,  1  <  j  <  4  can  be  expressed  in  the  form 

'i'-1  /  i  \i 

£  f— ^ 


j=i 


\2 ikj  aj 


(4.118) 


with  aj  :  R  — >  R,  1  <  j  <  7— 1  satisfying  condition  (4.110)  and  :  RxC+  — *  C 
satisfying  condition  (4.111).  Obviously,  I\  and  /4  are  already  in  the  form  (4.118). 
We  now  use  Lemma  3.3  to  show  that  /2,  h  can  also  be  expanded  in  the  form 
(4.118).  Observing  that  77(f)  =  0  for  all  t  £  (0,  Ti),  is  absolutely  continuous 
(see  Lemma  3.9),  and  that  af~i\  1  <  j  <  7  —  1  are  absolutely  continuous  (due 
to  the  assumption  of  the  induction),  we  can  use  formula  (3.5)  in  Lemma  3.3  to 
expand  /2  and  each  term  Js  (s  =  1, . . . ,  7  —  1)  of  /3  as 


h(t,k) 

Js(t,k) 


(4) 


(4.119) 

(4.120) 


dt(j-s-i)  MT)a(*-i) 


(2k)  (4121) 


with  bs  :  RxC+  — ►  C  uniformly  bounded  on  RxC+  (see  Lemma  3.3).  Therefore, 
J2  is  in  the  form  (4.118)  due  to  (4.119),  and  I3  is  of  the  form  (4.118)  due  to 
(4.121),  (4.116).  Thus,  mi+i(t,k)  can  indeed  be  written  in  the  form  (4.109) 
satisfying  conditions  (4.110),  (4.111).  □ 
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Corollary  4.11  Suppose  that  for  any  a>  0,  the  region  K(a)  C  C  is  defined  by 
the  formulae 

K(a)  =  {k\k€C,Im(k)>0,\k\>a}.  (4.122) 

Suppose  further  that  the  functions  :  R  x  C+  — *  C  are  defined  by  the 

formulae  (4-89),  (4-90),  (4-108)  and 


n^(t,k)  =  my(t,k)  —  — f  7?(r)e2!* ^  tW7(t,  (4.123) 

2tfc  Jo 

respectively.  Then  under  the  conditions  of  the  preceding  lemma,  there  exist  posi¬ 
tive  numbers  A,  ci,C2,c3  such  that 


\m(t,k)  -m^(t,k)\  < 

|n(*,fc)-n7(*,fc)|  <  |^ 


/or  all  ( t ,  k)  €  R  x  K(A),  and 

n{t,k)  c3 

m(t,k)  ~  |fc|'y 


(4.124) 

(4.125) 


(4.126) 


/or  a//  (<,&)  €  [Ti,oo)  x  K(A). 

Proof.  Due  to  (4.98),  the  norm  of  the  integral  operator  F*  in  (4.108)  is 
of  the  order  0(|&|-1)  for  any  k  €  C+.  from  which  we  observe  that  there  exists 
A  >  0,  such  that  (4.124)  is  true. 

Subtracting  (4.123)  from  (4.90),  we  obtain 

n(t,  k)  —  n7(f,fc) 

=  m(t , k)  —  m7(f,  k)  —  f  7/(r)e2l^<-r^(m(r,  k )  —  m7(r, k))di{.4.l21) 

ZZ  K  J 0 

Now,  the  estimate  (4.125)  is  a  direct  consequence  of  (4.127),  (4.124)  and  the 
fact  that  in  (4.127),  the  expression 

±p(r)t*k^  (4.128) 

is  uniformly  bounded  for  all  k  €  K(A),  — oo  <  r  <  t  <  oo. 

We  now  prove  (4.126)  by  showing  that  there  exists  a  positive  number  c3  such 
that 

n7(f,  k )  ^ 

m7(f,  fc) 


_C3_ 


(4.129) 
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for  all  (t,  k)  €  [7\,oo)  x  K(A).  According  to  Lemma  4.10,  m7(t,fc)  can  be 
expressed  in  the  form 


m7(t,  k)  =  1  +  £  (^)  “AO  +  (^)  <h(t,  k),  (4.130) 

with  aj,j  =  1,...,7  satisfying  conditions  (4.110),  (4.111).  Therefore,  we  cam 
assume  that  the  constant  A  has  been  chosen  such  that  for  all  (t,  k)  €  R  x  K(A), 

K»(*,*)l>£-  (4.131) 

Combining  (4.123)  with  (4.130),  we  obtain 

Tiy(t,  k)  -  m,y(t,  k)  +  I2(t ,  k)  +  h{t,  k)  +  I$(t,  k ),  (4.132) 

with  l2,h(t,k)  defined  by  (4.115),  (4.116),  and  h(t,k)  defined  by  the  formula 

/  1  \  *H-1  ft 

h(t,  k)  =  (^j)  jf  ,( tK(t,  k)e*W-’')jT.  (4.133) 

Noticing  that  tj(t)  =  0  for  all  t  >  T\,  we  have 

hlt.k)  =  (±)\«.k),  (4.134) 

=  (2- b,(t,k)  (4.135) 

for  all  (t,fc)  €  [Ti,  00)  x  K(A),  due  to  (4.119),  (4.121).  Consequently,  there 
exists  c  >  0  such  that 

I h(t,  k)  +  I3(t,  k)  +  Is(t ,  k)\  <  ^  (4.136) 

for  all  ( t,k )  €  pi,oo)  x  K(A ),  since  a-,(t,k),  bs(t,k)  are  bounded  for  all  (t,k)  € 
pi,  00)  x  K(A),  and  s  =  1, . . . , 7  —  1. 

Now,  (4.129)  follows  immediately  from  (4.132),  (4.136)  and  (4.131).  The 
estimate  (4.126)  is  a  direct  consequence  of  (4.129),  (4.124)  and  (4.125).  □ 

Lemma  4.12  Suppose  that  q  €  c2([0, 1]),  7  >  2,  is  absolutely  continuous 
and  q(x)  >  —1  for  all  x  €  R.  Then  there  exists  a  positive  number  c  such  that 

|p_(x,fc)-l|  <  jjL  (4.137) 


for  all  x  >  1,  k  £  C+ . 
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Proof.  According  to  Corollary  3.11  and  formula  (3.41), 

#-(*.*) 


p-{x,k)  = 


— ,  k ) 


(4.138) 


for  all  t  >  Ti  (i.e.,  for  all  x  >  1),  €  C+.  According  to  (4.89),  (4.90)  and 

(4.138) 

,  m(t,k) 

p-{x’k)=7(m  (4-i39> 

for  all  t  >  Ti,  €  C+.  Now,  the  lemma  follows  immediately  from  (4.139)  and 
(4.126).  □. 


Remark  4.13  By  a  similar  calculation,  one  can  show  that  under  the  conditions 
of  the  preceding  lemma,  there  exist  positive  numbers  A  >  0,  c  >  0  such  that 

-  1|  <  (4.140) 


for  all  x  <0,  k  €  C+ . 

Theorem  4.14  Suppose  that  q  €  Cq([0,1J),  q(x)  >  -1  for  all  x  £  R  and  q"  is 
absolutely  continuous.  Suppose  further  that 

D  =  {(x, k)\x  €  R,Im(k)  >  0}.  (4.141) 


Then 

(a)  (j>+  and  <j>-  are  continuous  functions  of  (x,  k)  and  analytic  functions  of  k  for 
all  x  €  R  and  k  €  C ; 

(b)  p+  and  p_  are  continuous  functions  of  (x,  k)  and  analytic  functions  of  k  in 

D; 

(c)  there  exists  a  positive  number  c  such  that  for  all  (x,  k)  6  D 

P+(x,k)  =  V1  +  -  i (T+'^T))  -  7k  +  *>■  (4-142) 

p-{x'k)  =  +  +  +  <4143> 
with  e+,e_  :  D  — *  C  continuous  functions  such  that 

\e+{x,k)\<-~,  (4.144) 

M*,*)l<j£jj-  (4.145) 
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Proof.  We  only  give  the  proof  for  <p_ ,  p_  since  the  proof  for  <f>+ ,  p+  is  identical. 
We  introduce  two  auxiliary  functions  <f>  and  <pi  via  the  formulae 


< fi(x,k )  =  4>_(x,k),  (4.146) 

<t>i(x,k)  =  (j>'_(x,k),  (4.147) 

so  that  the  equation  (2.1)  and  the  initial  condition  (2.5)  for  <f>-  can  be  rewritten 
as  a  system  of  linear  ODEs 

<t>'(x,k)  =  <j>x(x,k ),  (4.148) 

(f>[(x,k)  =  —  k2n2(x)<f>(x,k),  (4.149) 

subject  to  initial  conditions 

<f>{0,k)  =  1,  (4.150) 

4>i(0,k)  =  -ik.  (4.151) 


According  to  Lemma  3.6,  <j>,  (f>\  axe  continuous  functions  of  (x,  k)  and  entire 
functions  of  k  for  all  x  G  R  and  k  6  C,  from  which  (a)  follows  immediately. 
Similarly,  we  obtain  (b)  by  combining  (a)  with  (2.16)  and  the  fact  that  <£_(x,  k )  ^ 
0  for  all  (x,k)  €  D  (see  Remark  (4.3),  Lemma  4.4  and  Corollary  4.7). 

The  expansion  (4.143)  and  the  estimate  (4.145)  follow  immediately  from 
(3.41)  (see  Corollary  3.11  in  Chapter  3),  (4.83),  (4.84),  (4.87),  and  (4.88)  (see 
Lemma  4.9).  □ 

Corollary  4.15  Denote  by  p  p+  or  p_.  Then  under  the  conditions  of  the  pre¬ 
ceding  theorem,  there  exist  positive  number  Cj,C2  such  that 

e2ikftxp(r,k)dr  ^  (4152) 

for  all  t,x  €  [0, 1],  k  €  R,  or  for  all0<t<x<l,k€.  C+,  and 

\p'(x,k)\<c2,  (4.153) 

for  all  x  €  R,  k  €  C+ . 

Proof.  Due  to  Statements  (b),  (c)  of  Theorem  4.14,  the  real  part  of  the 
function 

2 ik  f  p(r ,  k)dr  (4.154) 

J  t 

is  uniformly  bounded  from  above  for  t,  x  €  [0, 1],  k  €  R,  or  for  all  0  <  t  <  x  <  1, 
k  €  C+,  from  which  (4.152)  follows  immediately.  Estimate  (4.153)  is  a  direct 
consequence  of  Statement  (c)  of  Theorem  4.14,  and  formulae  (3.50),  (3.51).  □ 

Global  upper  and  lower  bounds  for  the  impedance  functions  will  be  established 
in  Theorem  4.17.  We  first  obtain  a  partial  result  in  the  following  lemma. 
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Lemma  4.16  Suppose  that  for  any  positive  numbers  a,  a,  the  domain  K  (a,  a)  C 
C  is  defined  by  the  formula 

K(a,a)  =  {k\k  G  C,Re(k)  G  [— a,  a], Im(k)  G  [0, or]}-  (4.155) 

Then  under  the  conditions  of  the  preceding  theorem,  for  any  A  >  0,  there  exist 
positive  numbers  B,b,8  such  that 

|p+(x,A:)l  <  B,  (4.156) 

|p_(x,£)j  <  B ,  (4.157) 

Re{p+(x,  k))  >  b,  (4.158) 

i?e(p_(x,fc ))  >  b,  (4.159) 

in  the  domain  R  x  K(A,  8). 

Proof.  Since  the  proof  of  (4.156),  (4.158)  is  identical  to  that  of  (4.157), 
(4.159),  we  only  provide  the  latter.  Denoting  by  u ,  v  the  real  and  imaginary  parts 
of  p_  so  that 

p-(x,k)  =  u(x,k)  +  iv(x,k),  (4.160) 

the  Riccati  equation  (3.51)  for  p_  can  be  rewritten  in  the  form 

u'  =  -2  kuv,  (4.161) 

v'  =  —  k(v2  —  u2  +  n2),  (4.162) 

for  any  k  €  R.  Integrating  (4.161)  on  interval  [0,  x]  and  observing  that 

u(x,  k)  =  p_(x,  k)  =  1  (4.163) 

for  all  x  <  0,  k  G  C  (see  (3.53)),  we  have 

n(x,  k)  =  e~2k  ff  v(t’k)dt  >  o  (4.164) 

for  all  x,  k  G  R.  For  any  A  >  0,  p_,  u  =  f?e(p_)  axe  continuous  functions  of  (x,  k ) 
in  the  compact  domain  [0, 1]  x  K(A,8).  Therefore,  there  exist  positive  numbers 
bi,8,  B\  such  that 


u(x,k)  >  b\  >  0  (4.165) 

|p_(x,fc)i  <  Bu  (4.166) 

for  all  (x,  k)  G  [0, 1]  x  K (A,  8),  which  proves  the  estimates  (4.157),  (4.159). 

We  now  prove  the  estimates  (4.157),  (4.159)  for  all  x  >  1  using  the  formula 


1  —  bt(k)  +  i2b-(k)  s\n(kx  —  a_(fc )) 

1  +  bt(k)  +  2b-(k)  cos(kx  —  a_(fe))  ’ 


(4.167) 
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(see  Remark  (2.2)).  According  to  Remark  (2.2),  b(k)  >  0  is  a  real- valued 
continuous  function  of  k  €  C.  We  observe  that 

0  <  b{k)  <  1  (4.168) 

for  all  k  in  the  close  domain  K(A,8)  since  otherwise  if  b(k)  >  1,  the  real  part  of 
p-(l,k) 

u(\  M  =  _ 1  ~  _  (4.169) 

’  1  +  b2_(k)  -|-  26_(&)  cos(kx  —  at-(k)) 

will  be  non-positive,  contradicting  (4.165).  Due  to  (4.168),  (4.167),  there  exist 
positive  numbers  62?  B2  such  that 

u(x,k)  >  b2,  (4.170) 

\p.(x,k)\  <  B2 ,  (4.171) 

for  all  x  >  1,  k  €  K(A,8). 

Now,  (4.157),  (4.159)  follow  immediately  from  (4.165),  (4.166),  (4.170), 
(4.171),  and  (4.163).  □ 

Theorem  4.17  Suppose  that  q  €  Cq([0,  1]),  q(x)  >  —1  for  all  x  6  R  and  the 
second  derivative  of  q  is  absolutely  continuous.  Then  there  exist  real  numbers 
B  >  0,  b  >  0  such  that 


in  the  domain 


|p+(z,*)l 

< 

B, 

(4.172) 

|p_(x,*)| 

< 

B , 

(4.173) 

Re(p+(x,*)) 

> 

(4.174) 

Re{p-(x,k)) 

> 

b, 

(4.175) 

{(x,*)|x  €  R, 

Im(k )  >  0). 

(4.176) 

Proof.  Since  the  proof  of  (4.172),  (4.174)  is  identical  to  that  of  (4.173), 
(4.175),  we  only  provide  the  latter.  According  to  the  high-frequency  asymptotics 
(4.143)  in  Theorem  4.14,  there  exist  positive  numbers  A,  bi  such  that 

Re(p_{x,k))  >  (4.177) 

in  the  domain  D\  C  D  defined  by 

D\  =  {(x,  jfc)|x  €  R,  |*|  >  A,  Im(k)  >  0}.  (4.178) 

Since  p_(x,*)  is  a  continuous  function  of  (x,  k)  €  D\ ,  there  exists  a  positive 
number  B\  such  that 


|p_(x,*)|  <  Bu 


(4.179) 
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for  all  (x,fc)  €  D\.  For  such  a  number  A  >  0,  according  to  Lemma  4.16,  there 
exist  positive  numbers  8 ,  B2-,  b 2  such  that 

IM*,*)I  <  (4.180) 

Re(p-(x,k ))  >  b2,  (4.181) 

in  the  domain  D2  C  D  defined  by  the  formula 

D2  =  {(*,  *)|ar  €  R,Re(k)  €  [-A,A\Jm{k)  €  [0,*]}.  (4.182) 

Now,  according  to  Lemmas  4.5,  4.6,  there  exist  positive  numbers  B3,  b3  such 
that 

|p-(x,fc)|  <  B3,  (4.183) 

Re(p-(x,k))  >  63,  (4.184) 

in  the  domain  D3  C  D  defined  by 

D3  =  {(*,  k)\x  €  R,  Re(k )  €  [-A,  A],  Im(k)  >  <5}.  (4.185) 


The  estimates  (4.173),  (4.175)  for  (x,fc)  €  D  follow  immediately  from  the 

estimates  for  (x,  k)  €  Di,  D2,  D3  since  D  =  Z>i  U  D2U  D3.  □ 

The  following  theorem  furnishes  the  analytical  apparatus  for  the  error  analysis 
of  the  truncated  trace  formula  (see  (4.198)). 

Theorem  4.18  Suppose  thatq  €  CcP([0, 1]),  m  >  2,  is  absolutely  continuous 
and  q(x)  >  —1  for  all  x  G  R-  Then  there  exists  a  positive  number  a  such  that 

|p+(i,fc) -p_(x,fe)|  <  |~  (4.186) 


for  all  (x,  k)  €  R  x  C+ . 

Proof.  According  to  Lemma  4.12  and  Remark  (4.13),  (4.186)  is  true  for  all 
x  &  (0, 1).  In  order  to  prove  the  theorem  for  x  €  (0,1),  we  observe  that  pT  and 
p_  obey  the  same  Riccati  equation  (3.51)  due  to  (3.50),  (3.51).  The  difference, 
s  =  pT  —  p_,  satisfies  the  ODE 

s'(x,  k)  =  ik  (pT  +  p~)  s  (4.187) 


with  the  solution 


s(x,k)  =  e-ikfo(p+(t’k)+p-(t'k))dts(0,k). 


(4.188) 
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Corollary  4.15  indicates  that  there  exists  constant  b  >  0  such  that 

e-ikJ0x(JTiEk)+P-(t,k))dt  ^  (4.189) 

for  all  (x,  k)  €  [0, 1]  x  R.  Due  to  Remark  (4.13),  there  exists  a  positive  number 
c  such  that  for  all  k  G  R, 

|s(0,fc)|  =  \p+{0,k)  -p_(0,fc)|  =  |p+(0,fc)  -  1|  <  (4.190) 

Now,  (4.186)  for  x  €  (0, 1)  follows  immediately  from  (4.188),  (4.189),  (4.190). 

□ 


4.3  Trace  Formulae 

In  this  section,  we  prove  Theorem  4.19,  which  is  both  the  purpose  of  this  chapter, 
and  the  principal  analytical  tool  of  this  thesis.  Theorem  4.19  describes  the  so- 
called  trace  formulae  for  the  impedance  functions  p+,p~  (for  a  more  detailed 
discussion  of  the  term  ’’trace  formulae”,  see,  for  example,  [7]).  In  fact,  only  the 
formula  (4.194)  is  to  be  used  by  the  reconstruction  algorithm  of  the  following 
chapter.  We  present  the  formulae  (4.191),  (4.192),  (4.193)  for  completeness, 
since  some  of  them  appear  to  be  well-known,  and  attempts  have  been  made  to  use 
them  in  reconstruction  algorithms  (see,  for  example,  [8]).  See  also  Section  5.1 
below  for  a  more  detailed  discussion  of  the  use  of  trace  formulae  in  reconstruction 
schemes 


Theorem  4.19  (Trace  formulae)  Suppose  that  q  €  cJJ’QO,  1]),  m 
absolutely  continuous  and  q(x)  >  —  1  for  all  x  6  R.  Then 

>  2,  qW  is 

(a) 

1  fll 

\/l  +  ?(*)  =  Jjm,  —  J  ^  p+  (x,  k)dk. 

(4.191) 

(b) 

q'(x)  =  lim  —  (1  +  ?(x))  /  k-p+(x,k)dk. 
a— H-OO  la  J-a 

(4.192) 

(c) 

\A + ?(*)  =  ^  /a(p+(x> k ) + k))dk ■ 

(4.193) 

(d) 

cy  roo 

q'(x)  =  -(1  +  q(x))  I  (p+(x,  k)  -  p_(x,  k))dk. 

7T  J- oo 

(4.194) 
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More  precisely,  there  exist  positive  numbers  c1?  C2.C3.  c4  such  that 

yfi  +  q(x)  -  ^J_  P+  (x,  k)dk 
q'(x)  —  -r-(l  +  g(x))  f  k-p+(x,k)dk 

Id  J —a 

\Jl  +  q(x)  -  —  J_  (p+(x,  k)  +  p_(x,  k))dk 
q'(x)  -  -(1  +  q(x))  f  ( p+(x ,  k)  -  p_(x,  k))dk 

7T  J-a 

/or  a//  x  €  R. 


< 

£1 

5 

(4.195) 

a 

< 

£2 

(4.196) 

a 

< 

c3 

a2’ 

(4.197) 

< 

c4 

a(m_1)  ‘ 

(4.198) 

Proof.  Since  the  proofs  of  trace  formulae  (a),(b),(c),  and  (d)  are  similar,  we 
only  present  that  of  (d).  According  to  statement  (c)  of  Theorem  4.14,  there 
exists  c  >  0  such  that 


(p+(x,fc)  -p_(x,A:)) 


(  Si*)  1  \ 

V  2(1  +  <?(x))  ikj 


(4.199) 


for  all  (x,  k)  £  R  x  C+.  Denoting  by  T  the  upper  half  circle  of  radius  A ,  with 
clockwise  orientation,  in  the  complex  &-plane,  i.e., 


T  =  {*|*:  €  C+,  |*|  =  A},  (4.200) 

and  noting  that  p+  —  p_  is  an  analytic  function  of  k  G  C+,  we  obtain 

f  (p+(x,k)  —  p-(x,k))dk  =  f  (p+(x,k)  —  p-(x,k))dk.  (4.201) 

J—A  Jr 

Substituting  (4.199)  into  (4.201),  we  have 

J_A(p+{x,  k)  -  p-(x,  k))dk  =  ^l}q(xj)  +  °(k' ^  (4-202) 

from  which  (4.194)  follows  immediately. 

In  order  to  prove  the  estimate  (4.198),  we  rewrite  (4.194)  as 

q'(x)  =  -(1  +  q(x))  f  (p+(x,k)  -  p.(x,k))dk  +  1(a)  (4.203) 

7T  J-c 

with  1(a)  given  by  the  formula 

I(a)  =  “(!  +  ?(x))  +  Ja  )  (P+(x,k)  -  p.(x,k))dk.  (4.204) 
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Now,  formula  (4.74)  implies 

/(a)  =  -(1  +  q(x))  ^  +  Ja  (?+(*,  k)  -  p_(x,  k ))  dk , 

and  according  to  (4.186),  there  exists  a  constant  c4  such  that 

from  which  (4.198)  follows  immediately.  □ 


(4.205) 


(4.206) 


Chapter  5 

The  Reconstruction  Algorithm 


5.1  Reconstruction  via  trace  formulae-an  in¬ 
formal  description 

An  examination  of  the  formulae  (4.191)-  (4.194)  in  combination  with  the  Riccati 
equations  (3.50),  (3.51)  immediately  suggests  an  algorithm  for  the  reconstruc¬ 
tion  of  the  parameter  q  given  the  impedance  function  p+(x o,  k)  measured  at  some 
point  x0  6  R  outside  the  scatterer.  Namely,  one  is  tempted  to  substitute  one  of 
the  formulae  (4.191)-  (4.194)  (for  example,  (4.191))  into  (3.50),  obtaining 

P+{ x,  k)  =  -ik  (p2+{x,  k)  -  lim  ^  J  ^  p+{x,  k)dk^j  ^  ,  (5.1) 

and  attempt  to  view  (5.1)  as  a  differential  equation  for  the  function  p  :  Rl  x  R1  — ► 

C. 

Needless  to  say,  standard  existence  and  uniqueness  theorems  are  not  appli¬ 
cable  to  ‘differential  equations’  of  the  form  (5.1).  Furthermore,  in  order  to  be 
numerically  useful,  the  integral  in  (5.1)  would  have  to  be  replaced  with  some 
finite  quadrature  formulae.  The  latter  procedure  is  significantly  complicated  by 
the  fact  that  the  function  p+  is  defined  on  the  whole  real  line,  and  its  domain  of 
definition  has  to  be  truncated  before  discretization.  It  turns  out  that  the  solution 
of  (5.1)  is  not  unique,  except  in  a  very  carefully  chosen  class  of  functions  p.  Such 
a  class  of  functions  has  been  successfully  specified  (see,  for  example,  [8]).  The 
resulting  numerical  scheme  is,  however,  quite  expensive,  and  the  construction  is 
not  rigorous,  though  we  believe  that  this  could  be  made  so.  The  same  prob¬ 
lem  arises  if  one  attempts  to  use  the  trace  formulae  (4.192),  (4.193),  and  the 
conceptual  reason  for  this  situation  is  summarized  in  the  following  observation. 
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An  immediate  consequence  of  the  formula  (4.191)  is 


V1  +  «(*) 


lim  — 

a— *+oo  2 a 


p+(x,k)dk  +  [  p+(x, 
Jb 


(5.2) 


for  any  positive  real  b.  Thus,  the  'differential  equation’  (5.1)  can  be  replaced 
with 


p'+(x,  k)  =  —  ik 


lim  ~  l  [  p+(x,k)dk  +  j  p+(x,k)dk 
i— »+oo  2a  \J-a  Jb 


(5.3) 


and  a  convergence,  uniqueness,  etc.  proof  valid  for  (5.1)  would  also  be  valid  for 
(5.3),  unless  some  extremely  subtle  phenomenon  interfered. 

However,  given  a  smooth  scatterer  q,  for  any  e  >  0,  one  can  choose  a  suffi¬ 
ciently  large  b  that 


\Jl  +q(x)~P+{x,k) 


<  £ 


(5.4) 


for  any  k  >  b.  If  the  scattered  data  p+(x o,  k )  have  been  collected  at  some  point 
x0  outside  a  smooth  scatterer,  (5.4)  assumes  the  form 


|1  -p+(x,i)|  <  e. 


(5.5) 


In  other  words,  a  reconstruction  algorithm  using  the  ‘differential  equation’  (5.3) 
with  a  sufficiently  large  b  would  effectively  reconstruct  the  parameter  q(x)  for  all 
x  €  [0, 1]  from  a  single  measurement,  the  latter  being  equal  to  1  (!).  Another 
way  to  make  this  observation  is  to  notice  that  the  formula  (4.191)  is  simply  the 
WKB  approximation  to  the  impedance  function  p+ ,  and  that  in  the  WKB  regime, 
the  back-scattered  field  is  absent.  A  similar  problem  arises  if  one  attempts  to 
combine  formulae  (4.192),  (4.193)  with  (3.50),  and  view  the  result  as  a  “system 
of  ordinary  differential  equations” . 

In  the  case  of  a  discontinuous  scatterer  q,  the  WKB  expansions  (4.142), 
(4.143)  are  invalid.  On  the  other  hand,  the  trace  formulae  (4.191),  (4.194) 
are  valid  (if  the  limits  in  these  formulae  are  interpreted  properly),  and  can  be 
combined  with  the  equations  (3.50),  (3.51)  to  obtain  a  numerical  scheme  for 
detecting  discontinuities  in  the  scatterer.  If  q  is  piece-wise  constant,  such  a 
scheme  will  reconstruct  it  effectively,  and  time-domain  versions  of  this  procedure 
are  known  as  layer-stripping  algorithms  (see,  for  example,  [12],  [13],  [14]). 

While  the  author  failed  to  find  the  trace  formulae  (4.192),  (4.193)  in  the 
literature,  they  appear  to  be  well-known  among  specialists,  being  an  immediate 
consequence  of  the  WKB  analysis  of  the  equation  (3.50).  On  the  other  hand,  the 
formula  (4.194)  does  appear  to  be  new,  and  its  combination  with  the  equation 
(3.50)  immediately  leads  to  a  robust  reconstruction  algorithm.  While  we  post¬ 
pone  a  detailed  construction  and  analysis  of  such  a  scheme  till  Section  5.2,  in 
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the  following  observation  we  summarize  the  conceptual  reasons  for  its  analytical 
and  numerical  effectiveness. 

Formula  (4.198)  means  that  approximating  the  trace  formula  (4.194)  with 
its  ‘truncated’  version 

q'(x)  ~-(l  +  q(x))  [  {p+{x,k)-p.{x,k))dk,  (5.6) 

7 T  J —a 

we  make  an  error  of  the  order  where  m  is  the  smoothness  of  the  scatterer. 

Thus,  for  a  sufficiently  smooth  scatterer  and  a  sufficiently  large  a,  (5.6)  is  an 
extremely  good  approximation  to  the  trace  formula  (4.194). 

Now,  for  the  system  of  equations  (3.50),  (3.51),  (5.6),  it  is  not  hard  to 
prove  existence,  uniqueness,  etc.  theorems  of  the  type  valid  for  systems  of  ODEs 
(since  now  for  a  fixed  value  of  x,  the  functions  p+(x,  fc),p_(x,  k)  :  [—a,  a]  — »  C 
are  defined  on  a  compact  interval,  as  opposed  to  the  whole  line).  The  remainder 
of  this  thesis  is  devoted  largely  to  proving  such  facts  (see  Theorem  5.1  below), 
and  to  a  numerical  implementation  of  the  resulting  procedure.  The  latter  is  also 
quite  straightforward,  since  it  only  involves  constructing  a  quadrature  formula 
for  the  evaluation  of  the  integral  in  (4.194),  where  it  is  taken  over  an  interval 
of  finite  length.  Furthermore,  for  all  practical  purposes,  the  integrand  vanishes 
at  the  ends  of  the  domain  of  integration  together  with  all  its  derivatives,  com¬ 
pletely  obviating  the  issue  of  the  choice  of  the  quadrature  formula,  and  leading 
to  extremely  accurate  numerical  procedures  (see  Remark  6.3  below). 


5.2  Reconstruction  via  trace  formulae— a  for¬ 
mal  description 

Now,  we  are  prepared  to  construct  a  system  of  integro-differential  equations 
whose  initial  conditions  are  the  values  of  the  impedance  functions  p+,p~  mea¬ 
sured  outside  the  scatterer,  and  whose  solution  reconstructs  the  potential  q  for 


all  x  €  [0, 1).  We  will  consider  a  system  of  integro-differential  equations 

PU(*»*)  =  ~ik(P2a+(x’k)  -i1  +?«(*)))>  (5-7) 

Pa-{x,k)  =  ik(pl_(x,k)  -  (1  +  ?«(*))),  (5-8) 

?«(*)  =  7(1  +  ?«(*))/  (?«+(*,*) -p.-(*,*))<fe,  (5.9) 

TT  J— a 

with  respect  to  the  functions  pa+,pa_  :  [0,1]  x  [—a, a]  —>  C,  qa  :  [0,1]  — *  R , 
subject  to  the  initial  conditions 

Pa+(0,fc)  =  Po(fc),  (5.10) 

P«-(0,.*)  =  1,  (5.11) 

?(0)  =  0.  (5.12) 
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It  turns  out  that  for  sufficiently  large  a,  the  system  (5.7)-(5.12)has  a  unique 
solution  for  all  x  €  [0,1],  that  this  solution  is  stable  with  respect  to  small  per¬ 
turbations  of  the  initial  data  po(k)  ,  and  that  qa  converges  to  q  as  a  —*  oo.  The 
following  theorem  and  Lemmas  5. 2-5. 4  formalize  these  facts. 

Theorem  5.1  (Convergence  of  the  inversion  algorithm)  Suppose  that  q  €  c™([0, 1]) 
m  >  4,  q*'™')  is  absolutely  continuous  and  q(x)  >  —  1  for  all  x  6  R.  Then  there 
exist  constants  A  >  0,  c  >  0  such  that 

(q(x)  -  qa(x)l  < (5.13) 

for  all  x  €  [0, 1],  a  >  A. 

Since  the  proof  of  this  theorem  is  quite  involved,  we  break  its  technical  part 
into  three  lemmas  which  are  then  directly  used  in  the  proof  of  Theorem  5.1. 

Lemma  5.2  Suppose  that  q  €  ^([0, 1]),  m  >  4,  q ^  is  absolutely  continuous 


and  q(x)  >  —  1  for  all  x  €  R.  Suppose  further  that  the  function  space  E  is  defined 
by  the  formula 

E  =  {[a,/?,7]|a,/9  6  c([0,l]  x  [-a,  a]), 7  €  c([0, 1])},  (5.14) 

equipped  with  the  norm 

ll/H  =max(|H|,||^||,||7||),  (5.15) 

with  f  =  [a,/3,7]  €  E.  Finally,  suppose  that  for  any  a  >  0,  the  functions 
/a,u?,  ea  :  R  — >  R  are  defined  by  the  formulae 

fa(x)  =  -  /  {p+{x,k)  -p-(x,k))dk,  (5.16) 

7T  J -a 

w(x)  =  “(1 +  ?(*)),  (5-17) 

7T 

£a(x)  =  -w(x)  i^J  +J  ^  (p+(x,k)  -p-(x,k))dk.  (5.18) 

Then  the  error  function  u  =  [e+,e_,/i]  €  E  defined  by  the  formulae 

e+(x,Jc)  =  pa+(x,k)-p+(x,k),  (5.19) 

e_(x,&)  =  pa-(x,  k)  -  p_(x,  k),  (5.20) 

h(x)  =  qa(x)-q(x)  (5.21) 

satisfies  the  equation 

L{u)(x,k)  =  N(u)(t,k)  +  [0,0,  £.(<)],  (5.22) 
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where  L,  N  :  £  —*  £  are  defined  by  the  formulae 

e+(x,  k)  -  ikf0x  h(t)e-2'kf'Zp+{T’k]dTdt 
L(u)  =  e_(x,k)  +  ikf0xh(t)e2ikf'Xp-(T’k)dTdt  (5.23) 

.  h(x)  -  fo  (h{t)fa(t)  +  w(t)f“a(e+  -  e~)(t,  z)dz^j  dt  _ 

-ikfZel(t,k)e-2ikftXp+(r’k)dTdt 

N(u)  =  ik  f*  e2_(t,k)e2tkf*  p~{T'k)dTdt  •  (5-24) 

.  f /o*M*)/-«(e+(M)  -  e_(t,  z)dzdt  _ 

Proof.  We  know  that  the  functions  p+,p~,q  satisfy  the  ODEs 

p'+(x ,  A;)  =  -ik(p2+(x,  k)  -  (1  +  9(x)),  (5.25) 

p'_(x,k)  =  ik(p2_(x,  k)  -  (1  +  g(x)),  (5.26) 

q\x)  =  -(l+q(x))  f  (p+(x,k)  -P-(x,k))dk,  (5.27) 

7 r  j-oo 

for  all  x€R,  any  k  €  C+,  and  that  the  functions  pa+,pa-,qa  satisfy  the  ODEs 

p'a+(x,k)  =  -ik(p2a+(x,k)-(l  +  qa(x))),  (5.28) 

P«-(Xik)  =  ik(p2a_(x,k)  -  (1  +  90(x))),  (5.29) 

q'a(x)  =  f(!  +  <l«(x))  J_a(P*+(x’z) -pa-[x,z))dz.  (5.30) 

for  all  (x,  fc)  €  [0, 1]  x  [—a,  a],  subject  to  initial  conditions 

jPa+(0i  k)  =  p+(0,k)i  (5.31) 

pa-(0,  k)  =  p_(0,fc)  =  l,  (5.32) 

9.(0)  =  ?(0)  =  0.  (5.33) 

for  all  k  €  [—a,  a].  Subtracting  equations  (5.25),  (5.26),  (5.27)  from  equations 
(5.28),  (5.29),  (5.30)  respectively,  we  observe  that  [e+,  e_,  h]  (see  (5.19),  (5.20), 
(5.21))  satisfies  the  ODEs 

e+(x,fc)  =  —  ik  (2j?+(x,  k)e+(x,  k)  +  e+(x,  k)  —  h(x)')  ,  (5.34) 

e'_(x,fc)  =  ik  (2p-(x,  fc)e_(x,  k)  +  ei(x,  k)  —  h(xfj  ,  (5.35) 

h’(x )  =  h(x)fa(x)  +  w(x)  J  (e+(x,z)  —  e_(x,x))dz 
2  ra 

4 — h(x)  /  (e+(x,  z)  —  e_(x,  z))dz  +  ea(x),  (5.36) 

7T  J— a 

subject  to  the  initial  conditions 


e+(0,  k)  =  e_(0,  k )  =  h(  0)  =  0. 


(5.37) 
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We  now  convert  the  initial  value  problem  (5.28)-(5.33)  as  a  system  of  integral 
equations.  Multiplying  (5.34)  by  the  function 

e2ikf*p+(t,k)dt'  (53g) 


we  have 

A  (e2ikfo  r+VWe+ix,  k ))  =  -ik  ■  e2ikfo  P+Mdt  (e*  (x,  k)  -  h(x))  .  (5.39) 

Integrating  the  result  over  the  interval  [0,  x],  we  obtain 

e+(x,  k)  —  ik  fXh(t)e-2ikf>+{T'k)dTdt  =  -ik  f  e\ (t,  k)e~2ikS' p+{T'k)drdt.  (5.40) 
Jo  Jo 

A  similar  calculation  reduces  (5.35)  to  the  equation 

e_(x,  k)  +  ik  f  h{t)e2ikl •  p~^)drdt  =  ik  fX  g2  kyikft*  P-(r,k)drdt  (5,41) 
Jo  Jo 

An  integration  of  (5.36)  over  the  interval  [0,  x]  converts  (5.36)  into  the  integral 
equation 

h{x)  -  jf  h(t)fa(t)dt  -  £w{t)L  (e+(t,z)  —  e-(t,z))dzdt 

=  —  /  h{t)  f  (e+(t,z)  —  e-(t,  z)dzdt  +  /  ea(t)dt.  (5.42) 

7 T  Jo  J-a  Jo 

Clearly,  equations  (5.40),  (5.41),  (5.42)  is  equivalent  to  (5.22),  which  completes 
the  proof.  □ 

Lemma  5.3  Under  the  conditions  of  Lemma  5.2,  there  exists  a  positive  number 
ci  such  that  for  any  f,g  €  S,  there  exist  continuous  functions  61,82  :  [0,1]  x 
[—a,  a]  — *  C,  83  :  [0, 1]  — *■  C  such  that 

N(f)(x,  k )  -  N(g)(x,  k )  =  [<5i(x,  k),  h2(x,  k),  [  83(t)dt],  (5.43) 

Jo 

and 

max(M,  ||62||,  ll^sl!)  <  ci  •  a  •  max  (||/||,  \\g\\)  \\f  -  g\\.  (5.44) 

Proof.  Formula  (5.43)  is  a  direct  consequence  of  (5.24).  In  fact,  we  have 

81(x,k)  =  -ik  [X(f2(t,k)-g2(t,k))e-2ikI>^dTdt,  (5.45) 

Jo 

62 (x,k)  =  ik  (5.46) 

Jo 

Hx)  =  “  /  {fsWifi- f2)(x,z)-g3(x){g1-g2){x,z)}dz  (5.47) 

7T  J-a 
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for  any  /  =  /3]  €  E,  g  =  [gi,  92,93]  €  E.  In  order  to  prove  (5.44),  we  first 

observe  that  due  to  Corollary  4.15,  there  exists  a  positive  number  c4  such  that 


e~2ik  ft  p+(T’k)dr 


<  c4, 


(5.48) 


e2ik  ft  p-lT’k)dT 


<  c4, 


(5.49) 


for  all  t,x  6  [0,1],  k  G  R.  Observing  that  \k\  <  a,  0  <  x  <  1,  and  using  the 
estimate  (5.48),  we  obtain  the  estimate 


INI  <  a  •  <*11/  +  g ||  •  11/  -  9 1|-  (5.50) 

A  similar  calculation  shows  that 

INI  <  a  ■  c4\\f  +  g\\  •  ||/  -  $||,  (5.51) 

and  we  obtain  the  estimate  for  63  by  first  regrouping  (5.47): 

INI  =  -  sup  |((/3(*)  -g3{t))  f  (/l(M)-  f2(t,z))dz 
*€[0,1]  IV  J-a 

+«3(i) 

<  l(\\f~9\\  Jj\f  +  9\\dz+y\\  J°nf~9\\dz) 

<  «- -I1/-ffll(ll/  +  sil  +  2||ff||)  (5.52) 

7T 

Now,  (5.44)  follows  immediately  from  (5.50),  (5.51),  (5.52).  □ 


Lemma  5.4  Under  the  conditions  of  Lemma  5.2,  there  exist  positive  numbers 
C2,  c3  such  that  for  any  6  <E  E  of  the  form 


6(x,  k )  =  Nx,  k ),  82(x,  k),  f  S3(t)dt], 

Jo 

(5.53) 

the  linear  equation 

II 

(5.54) 

has  a  unique  solution  v  =  [i>i,V2,u3]  £  £.  Furthermore, 

H  <  c2  •amax(||<51||,P2|j)  +  c3||<53||. 

(5.55) 

Proof.  We  only  need  to  prove  (5.55),  since  the  existence  and  uniqueness  of 
the  solution  v  of  the  linear  equation  (5.54)  is  a  direct  consequence  of  the  estimate 
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(5.55).  Due  to  (5.23),  (5.53),  the  equation  (5.54)  can  be  rewritten  in  the  form 


vx{x,k)  =  ik  fXv3{t)e-2,k^p^T'k)dTdt  +  81(x,k), 

Jo 

v2(x,  k)  =  -ik  f  v3(t)e2ikf*p-(r'k)dTdt  +  S2(x,  k) 

Jo 

^(x)  =  V3  (t)fa(t)dt  +  Sow{t]L  z )  —  v2(t,  z)dzdt 

+  /  S3(t)dt. 

Jo 


(5.56) 

(5.57) 


(5.58) 


In  order  to  prove  (5.55),  we  first  eliminate  t>i ,  v2  from  (5.58)  and  obtain  an 
estimate  for  v3.  Subtracting  (5.57)  from  (5.56),  and  integrating  the  result  over 
the  interval  [—a,  a],  we  obtain 

f  (vi(x,z)  -  v2(x,z))dz 

J —a 

=  £  V3{t)  ja  iz  +  e2 izftxP-(r,z)dr"j  dzdt 

+  [  {di(x,z)-82{x,z))dz 
J  —a 

=  f  ga(xj)v3(t)dt  +  4a- $a(x),  (5.59) 

Jo 


with  ga  :  [0, 1]  x  [0, 1]  —*  C,  sa:  [0, 1]  — ►  C  given  by  the  formulae 


ga(x,  i)  =  f  iz  (W*  p+(T’*)dT  +  e2,zf'  p"(T-*)dT)  dz, 
1  fa 

Sa (*)  =  a(8i{x,  z )  -  62(x,  z))dz. 


(5.60) 

(5.61) 


Combining  (5.58)  with  (5.59),  we  obtain 


U3(z)  =  J  v3(t)fa(t)dt  +  j  w(t)  J  ga(t,r)v3(r)drdt 
+4a  J  w(t)sa(t)dt  +  J  63(t)dt. 


(5.62) 


We  will  obtain  the  estimate  (5.72)  for  v3  (see  below)  by  first  proving  (5.63), 
(5.64),  (5.65),  and  (5.69)  for  functions  fa,w,ga,sa.  Obviously,  there  exist  con¬ 
stants  c5  >  0,  C6  >  0  such  that 


for  all  x  €  R  due  to  (5.17),  and 


Mx)|  <  c5 


!/■(*)!  <  Cg 


(5.63) 


(5.64) 
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for  all  x  £  [0, 1],  any  a  >  0,  due  to  (5.16),  (4.198),  and 

M*)|  <max(||6i||,||62||).  (5.65) 

due  to  (5.61).  Observing  that 

£  2  •  e~2i2f<X =-fz-  e2izf*t^dTdz  (5.66) 

due  to  (4.72),  and  combining  (5.60)  with  (5.66),  we  have 

ga{x,t)  =  Ja  iz  ^e2iz  ft  p~U'z)dT  _  e2»*  J7  P+(T,z)dr^  dz  (5  67) 

According  to  Theorem  4.18,  for  any  x  e  R,  the  function 

P-(x,k)  -p+(x,k)  (5.68) 

decays  uniformly  like  k~m ,  for  k  £  R,  and  consequently,  the  integrand  in  (5.67) 
decays  like  /r(m~2)  uniformly  with  respect  to  t,x  £  [0, 1].  Since  we  have  assumed 
that  m  >  4,  there  exists  a  constant  c7  >  0  such  that 

<c7,  (5.69) 

for  all  t,x  £  [0, 1],  a  >  0.  Now,  combining  the  integral  equation  (5.62)  with  the 
estimates  (5.63),  (5.64),  (5.65),  (5.69),  we  have 

b3(x)l  <  Cg  \v3(t)\dt  +  c5  •  c7  J  J  \v3(r)\dtdr 
+4a  •  c5  •  x max (||$i ||,  ||£2||)  +  INI 
-  [  (C6  +  c5c7)b3(0l^  +  4a-c5xmax(!|(51||,||52||)  +  ||d3||  (5.70) 

J  0 

Now,  the  estimate  for  v3  follows  from  Gronwall’s  inequality  (see  Lemma  3.5), 
M*)l  <  4a-c5-xmax(p1||,||62||)+  ||<53|| 

+  c&jo  (4a  •  c5t max  (ll^U,  ||62jj)  +  ||^3j|)  e^x-<^C6+CsC7>dt,  (5.71) 

with  c8  =  C6  +  c5  •  c7.  Clearly,  there  exist  positive  numbers  eg,  Ci0  such  that 

M*)l  ^  eg  max  (li^  ||,  ||<52||)  +  doll^aJI,  (5.72) 

for  all  x  £  [0, 1],  we  thus  have  the  estimate  for  t>3  (see  (5.55)). 

In  order  to  obtain  similar  estimates  for  ua,  v2,  we  first  provide  an  estimate  for 
the  derivative  of  v3.  Differentiating  (5.62),  we  have 

t4(x)  =  v3{x)fa{x)dt  +  xv(x)  [  ga(x,  t)v3{t)dt  +  4o  •  tx;(x)sQ(x)  +  63(x).  (5.73) 

J  0 
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Combining  (5.73)  with  (5.72),  we  observe  that  there  exist  positive  numbers 
cn,ci2  such  that 


KC*)!  <cn  max  (II^IIJ^ID  +  ciiPsli.  (5.74) 


Integrating  by  parts  in  (5.56)  yields 


■*<*•*>  - 

V3(x) 


=  6i(x,k)+- 


2  \p+(x,k) 

_  rx  ^(Qp+fo  fc)  -  «s k)__2ik  [;P+(T,k)dTdi 

Jo  p\  (t,  k)  j 


.  (5.75) 


For  all  (x,  k)  €  R  x  C+,  p+  is  uniformly  bounded,  Re(p+)  is  uniformly  bounded 
from  below  by  a  positive  number  (see  Theorem  4.17).  Due  to  Corollary  4.15, 
p'+  and 

e-2ik  ftx  P+(r,k)dr  (5  ?6) 

are  uniformly  bounded  for  all  x,f  €  [0,1],  k  €  R-  Therefore,  combining  (5.75) 
with  (5.72),  (5.74),  (4.152),  we  observe  that  there  exist  positive  numbers  ci3,  ci4 
such  that 

Mz,&)l  <  ci3 max (||6i||,  H^ll)  +  Ci4||63||j  (5.77) 

and  a  similar  calculation  shows  that 


\v2(x,k)\  <  C13 max(Pi||,  p2||)  +  1|<53||  (5.78) 

for  all  x  €  [0, 1],  k  €  [—a,  a].  Now,  the  estimate  (5.55)  follows  immediately  from 
(5.72),  (5.77),  (5.78).  □ 

Using  Lemmas  5.2,  5.3,  5.4,  we  now  proceed  with  the  proof  of  Theorem 
5.1. 

Proof  of  Theorem  5.1.  Theorem  4.19  implies  that  there  exists  positive 
numbers  61,62  such  that 


\fa(x)\  <  6„  (5.79) 

M*)l  S  <5-8°) 

We  prove  the  theorem  by  showing  that  there  exist  positive  numbers  A ,  c  such 
that  for  all  a  >  A,  the  solution  u  =  [e+,  e_,  h]  €  £  exists  (see  Lemma  5.2  for  the 
definitions  of  u,E),  and  that 
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We  will  obtain  the  solution  u  of  the  equation  (5.22)  via  the  following  iterative 
procedure: 

uq  =  0,  (5.82) 

L(un+ 1)  =  N(un)  +  [0,0,  ea(t)],  (5.83) 

with  L ,  N  defined  by  the  formulae  (5.23),  (5.24),  respectively.  Clearly,  we  only 
need  to  show  that  there  exist  positive  numbers  A,  c  such  that  for  all  a  >  A 

IWI  S  (5.84) 

and  the  sequence  un,  n  =  0, 1, . . .  converges  (to  the  solution  u).  The  first  iterate 
Ui  satisfies  the  equation 

£(Ul)  =  [0,0,  ea(<)],  (5.85) 

and  according  to  (5.18),  (4.198),  there  exists  a  constant  b  such  that 

IWI  ^  JrTy  (5-86) 

Combining  (5.85)  with  (5.55)  and  (5.86),  we  observe  that  there  exist  constant 
c4  such  that 

IM  <  (5-87) 

Now,  we  choose  a  constant  A  >  0  such  that 

a  •  Ci(c2  ■  a  +  c3)||ui||  <  -  (5.88) 

for  all  a  >  A.  Defining  «_i  =  2u\  for  convenience,  we  prove  by  induction  that 

|K+l~Un!|  <  -  Un-iH,  (5.89) 

IK+ill  <  2||tt1||,  (5.90) 

for  all  n  >  0,  a  >  A. 

The  case  n  =  0  is  a  trivial  one.  For  n  >  1,  (5.83)  indicates  that 

L(un+i  -  un)  =  N(un)  -  N(un-i),  (5.91) 

Due  to  Lemma  5.3,  there  exist  continuous  functions  61,62  ■  [0, 1]  x  [—a,  a]  — ♦  C, 
63  :  [0, 1]  — *  C  such  that 

N(un)  —  N(un~i)  =  [<5i(x,fc),<52(x,fc),  /  63(t)dt], 


(5.92) 
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Now,  combining  (5.91),  (5.92)  with  (5.55),  (5.44),  and  the  assumption  of  the 
induction,  we  obtain 

||«n+i-wn||  <  C2-amax(||6i||,||£2||)  +  c3|NI 

<  a  ■  ci(c2  •  a  +  c3)max(||un||,  ||un-i||)  ||u»  ~  «»-i|| 

<  i||tin  (5.93) 

which  proves  (5.89).  The  estimate  (5.90)  is  a  direct  consequence  of  (5.89). 

Finally,  the  sequence  un,  n  =  0,1,...  converges  to  the  solution  u  due  to 
(5.89),  and  therefore 

INI  S  (5.94) 

for  all  a  >  A  due  to  (5.90),  (5.87),  which  was  to  be  proved.  □ 

Remark  5.5  The  proof  above  requires  that  q  €  cm(R),  with  m  >  4.  At  the  ex¬ 
pense  of  a  considerable  increase  in  the  complexity  of  the  proof,  it  is  not  difficult  to 
extend  this  result  to  m  >  2.  However,  our  numerical  experiments  (see  the  follow¬ 
ing  chapter)  indicate  that  the  scheme  works  quite  well  for  continuous,  piecewise 
continuously  differentiable  q,  and  even  for  piecewise  continuously  differentiable  q 
with  finite  number  of  jumps.  In  the  latter  two  cases,  the  rates  of  convergence  of 
the  algorithm  are  1/a  and  l/y/a,  respectively. 


Chapter  6 


Implementation  and  Numerical 
Results 


6.1  Implementation 

In  implementing  the  algorithm  of  this  thesis  (see  Chapter  5.2),  the  integral 

f  (p+(x,Jc)-p_(x,k))dk  (6.1) 

J —a 

in  equation  (5.9)  is  approximated  by  the  trapezoidal  sum 

M- 1 

T{h)  =  h  (P+(x,kj) ~P-ix,kj)) 

j=—M+l 

+^((p+{x,  -°)  -  P-(x,  -a))  +  (?+(*>  a)  ~  P-{x,  a))),  (6.2) 

with  h  =  a/M,  k:  —  jh ,  j  =  — M , . . . ,  M.  Since  for  real  k,  p+(x,  —k)  =  p+(x,  k), 
P-(x,—k)  =  p_(x,&)  (see  Observation  (4.8)),  the  ODEs  (5.7),  (5.8),  (5.9)  are 
discretized  in  the  fc-space  using  M  +  1  nodes  kj  =  jh,  j  —  0, . . . ,  M,  leading  to 
a  system  of  2 M  +  3  ODEs 

ki)  =  ~ik j  (p/U(z>  k i)  -  (i  +  Qh(x)))  ,  (6.3) 

kj)  =  *'*j  (pL(*>  *j)  -  (1  +  qh{x)j)  ,  (6.4) 

Ah  (M~  i 

?*(*)  =  —  (l  +  9fc(*))  ( 

+  ^i?e{Pft+(M)  -pft_(x,0)  +  pfc+(x,a)-pft_(x,a)})  (6.5) 
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subject  to  the  initial  conditions 


Ph+(  0,kj) 

=  Po(kj), 

(6.6) 

Ph-iO.kj) 

=  1, 

(6.7) 

qh(  o) 

=  0 

(6.8) 

(see  (5.10)-  (5.12)).  These  ODEs  are  then  solved  using  a  standard  4-th  order 
Runge-Kutta  scheme. 

When  an  integral  is  discretized  via  a  quadrature  formula,  the  rate  of  conver¬ 
gence  of  the  quadrature  is  critical  to  the  numerical  performance  of  the  algorithm. 
It  turns  out  that  while  the  estimate 


p+(x,k) -p_(x,k)  =  0(a  m)  (6.9) 

(see  Theorem  4.18)  ensures  a  rapid  convergence  of  qa  to  q  as  a  grows  (see  Theorem 
5.1),  it  also  guarantees  a  rapid  convergence  of  the  trapezoidal  quadrature  (6.2) 
to  the  integral  (6.1).  This  fact  is  formalized  in  the  following  lemma.  Its  proof  is 
based  on  the  Euler-Maclaurin  summation  formula  (see,  for  example,  [10]),  and  is 
omitted,  since  it  is  quite  involved,  and  incidental  to  the  purpose  of  this  thesis. 


Lemma  6.1  Suppose  that  q  €  c™([0, 1]),  m  >  2,  <^m)  is  absolutely  continuous 
and  q(x)  >  —1  for  all  x  €  R.  Then  there  exist  positive  numbers  c„,  n  =  0, ... 
such  that 


dn  (p+(x,a)  -p_(x,a)) 


dkn 


< 


w 


(6.10) 


Furthermore,  for  any  f3  >  0,  b  >  0,  there  exists  a  constant  c  >  0  such  that 


f  (p+(x,k)-p-(x,k))dk-T 

J — a 


(6.11) 


Using  the  estimate  (6.11),  and  reproducing  the  proof  of  Theorem  5.1  almost 
verbatim,  one  can  prove  the  following  theorem. 

Theorem  6.2  Suppose  that  q  €  c™([0, 1]),  m  >  4,  q ^  is  absolutely  continuous 
and  q(x)  >  —1  for  all  x  €  R.  Suppose  further  that  for  given  r  >  0,  s  >  0. 
q(r,s,x )  denotes  the  solution  qh  of  the  system  (6.3)-  (6.8)  with  h  =  r/as.  Then 
for  any  a  >  0,  /3  >  0,  there  exist  constants  A  >  0,  c  >  0  such  that 

\q{x)  -  q(a,  /?,  z)|  <  (6.12) 


for  all  x  €  [0, 1],  a  >  A 


50 


CHAPTER  6.  IMPLEMENTATION  AND  NUMERICAL  RESULTS 


6.2  Numerical  Results 

We  have  applied  the  algorithm  of  this  thesis  to  the  reconstruction  of  several  types 
of  scatterers,  from  infinitely  differentiable  q  to  discontinuous  q.  The  computations 
were  performed  in  double  precision  on  a  SPARC  1  computer  without  the  use  of 
the  accelerator.  The  results  of  numerical  experiments  for  four  classes  of  scatterers 
are  presented  in  this  chapter. 

In  the  first  class  (Examples  1-2.2)  are  scatterers  satisfying  the  smoothness 
conditions  of  Theorem  6.2.  In  the  second  class  (Example  3)  is  a  scatterer  q  vio¬ 
lating  the  smoothness  conditions  only  mildly  (it  is  continuous,  but  its  derivative 
is  discontinuous  at  two  points).  In  the  third  class  (Examples  4.1,  4.2)  are  scat¬ 
terers  that  strongly  violate  the  smoothness  conditions  by  being  discontinuous. 
Finally,  in  Example  6  a  scatterer  with  an  index  of  diffraction  that  changes  in 
several  order  of  magnitude  is  reconstructed.  As  is  well-known,  scatterers  of  this 
type  are  difficult  to  recover  due  to  strong  back  scattering. 

We  also  performed  a  crude  test,  in  Example  5,  of  stability  of  the  algorithm 
by  truncating  the  scattering  data  po(kj ),  j  =  1, . . . ,  M  after  1, 2,  or  3  digits.  The 
truncated  scattering  data  are  subsequently  used  in  reconstructions. 

In  Tables  6.1-  6.9,  hk  denotes  step  size  of  the  trapezoidal  rule  in  the  ^-interval 
[0,a],  Nx  denotes  the  number  of  points  in  the  x-interval  [0,  2tt],  E2,E°°  repre¬ 
sent  the  relative  L2  and  maximum  norm  of  error  of  the  reconstructed  scatterer, 
respectively.  In  Figures  6.1-  6.15,  dotted  lines  denote  the  exact  solution,  while 
solid  lines  denote  the  numerical  reconstruction.  In  all  examples,  for  a  given  a,  hk 
and  Nx  were  chosen  such  that  further  decrease  of  hk  and  increase  of  Nx  brought 
no  improvements  on  the  accuracy  of  the  reconstruction. 

Remark  6.3  In  order  to  obtain  the  scattering  data  p+( 0,k)  for  the  Examples 
1-3,  the  scattered  field  4>Scat+  was  obtained  as  a  solution  of  the  boundary  value 
problem  (2.10),  (2.8),  (2.9)  via  a  high  order  algorithm  described  in  [15].  The 
parameters  in  the  scheme  were  chosen  in  such  a  manner  that  at  least  14-digit 
accuracy  was  always  maintained.  Formulae  (2. 4),  (2.30)  were  then  used  to 

obtain  p+(0,  k )  from  4>SCat+- 

In  Examples  4-1  and  4-2,  a  standard  procedure  for  the  solution  of  the  initial 
value  problem  (2.1),  (2.26)  (for  <p+ )  with  piecewise  constant  q  was  used  (see,  for 
example,  [16]).  Here,  the  solutions  were  obtained  with  at  least  15  correct  digits. 
The  scattering  data  p+(Q,k)  were  obtained  from  (p+  via  formula  (2.30). 

Remark  6.4  In  the  examples  below,  no  effort  was  made  to  optimize  the  code 
used,  either  from  the  algorithmic  or  from  the  programming  point  of  view.  For 
example,  we  used  the  Runge-Kutta  scheme  to  solve  ODEs  (6.3),  (6-4),  (6.5). 
While  it  produced  satisfactory  results  in  our  experiments,  it  is  by  no  means  the 
most  efficient  scheme  for  the  solution  of  problems  of  this  type. 
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a 

hk 

Nx 

E2 

£OC 

t  (sec.) 

5 

0.1 

80 

0.146  xlO-2 

0.153  xlO"'2 

0.600 

10 

0.1 

300 

0.354  xlO-5 

0.415  xlO-5 

4.41 

10 

0.05 

600 

0.177  xlO"5 

0.183  xl0~5 

16.7 

10 

0.05 

1200 

0.175  xlO"5 

0.184  xlO-5 

34.2 

20 

0.05 

2400 

0.759  xlO'9 

0.108  xlO"8 

141 

20 

0.05 

4000 

0.988  xlO"10 

0.143  xlO-9 

235 

20 

0.025 

4000 

0.982  xlO"10 

0.142  xlO"9 

498 

Table  6.1:  CPU  Times  and  Accuracies  for  Example  1 


Example  1.  Reconstruction  of  a  Gaussian  distribution 

q(x)  =  e~(~ ) 

where  the  variant  a  given  by  the  formula 

a  =  =  0.5175854235 . . . 


(6.13) 


(6.14) 


was  chosen  such  that  the  function  is  effectively  zero  to  double  precision  outside 
the  interval  [0,  2tt).  The  results  of  this  numerical  experiment  are  depicted  in 
Table  6.1  and  Figure  6.1.  For  all  practical  purposes,  the  scatterer  (6.13)  is 
a  c°° -function  in  R  with  the  support  on  the  interval  [0, 27r],  and  therefore  the 
algorithm  converges  extremely  rapidly,  as  demonstrated  in  Figure  6.1  where  the 
two  graphs  of  the  exact  q  and  the  reconstructed  q  are  almost  identical. 


In  the  following  two  examples,  we  reconstruct  oscillatory  scatterers  of  the 
form 


?(*)  =  2  ci(i  “  cos(njx)), 

i=i 


(6.15) 


with  nj,Cj,  j  —  1,2,3  given  below.  For  given  the  coefficients  c:  were  chosen 
in  such  a  manner  that  q  is  five  times  continuously  differentiable  for  all  x  €  R, 
so  that  the  rapid  convergence  of  the  reconstruction  algorithm  is  guaranteed  (see 
Theorems  5.1,  6.2) 


Example  2.1.  A  less  complicated  scatterer  is  given  by  the  formula 


q{x)  =  0.3  ((1  -  cos(2x))  -  ^(1  -  cos(3x))  -I-  ^(1  -  cos(4x)))  .  (6.16) 


Reconstructions  were  performed  with  a  =  7,  14.  The  results  of  this  experiment 
are  depicted  in  Table  6.2  and  Figure  6.2  Since  the  scatterer  is  smooth,  p+(x,  k)  — 
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Table  6.2:  CPU  Times  and  Accuracies  for  Example  2.1 
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Figure  6.2:  Reconstruction  of  Example  2.1  with  a  =  7 


Figure  6.3:  Real  Part  of  the  Scattering  Data  po(k)  in  Example  2.1 
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a 

hk 

Nx 

E2 

E°° 

t  (sec.) 

10 

0.1 

300 

0.288  xlO"1 

0.376  xlO'1 

4.41 

10 

0.025 

600 

0.281  xlO-1 

0.367  xlO'1 

35.3 

10 

0.025 

1200 

0.281  xlO"1 

0.367  xlO-1 

70.4 

20 

0.1 

400 

0.395  xlO"2 

0.754  xlO-2 

11.4 

20 

0.025 

800 

0.127  xlO"2 

0.226  xlO-2 

98.7 

20 

0.025 

1600 

0.127  xlO-2 

0.220  xlO"2 

197 

40 

0.025 

800 

0.788  xlO"4 

0.300  xlO"3 

202 

40 

0.025 

1600 

0.878  xlO"5 

0.290  xlO"4 

404 

Table  6.3:  CPU  Times  and  Accuracies  for  Example  2.2 


a 

Nx 

E2 

E°° 

t  (sec.) 

5 

0.1 

75 

0.482  xlO-1 

0.829  xlO'1 

0.590 

10 

0.1 

150 

0.239  xlO'1 

0.462  xlO'1 

2.19 

20 

0.1 

300 

0.119  xlO'1 

0.283  xlO'1 

8.47 

Table  6.4:  CPU  Times  and  Accuracies  for  Example  3 


p_(x,£)  decays  rapidly  as  k  grows.  In  particular,  the  scattering  data  Re(po(k)) 
approaches  1  rapidly,  as  can  be  seen  in  Figure  6.3. 

Example  2.2.  A  more  complicated  scatterer  is  given  by  the  formula 

q(x)  =  0.4  ((1  -  cos(3x))  -  p^(l  ~  cos(llx))  +  ^(1  -  cos(12x))  .  (6.17) 

Reconstructions  were  performed  with  a  =  10,  20.  The  results  of  this  experiment 
are  depicted  in  Table  6.3  and  Figure  6.4. 

Example  3.  In  this  example,  we  reconstruct  a  scatterer  defined  by  the 
formula 

^(x)  =  0.2  •  sin(x).  (6.18) 

Note  that  q 1  is  discontinuous  at  the  points  x  =  0,  27t,  and  as  a  result  p+(x,  k )  — 
p_(x,fc)  decays  like  1/A:,  as  can  be  seen  in  Figure  6.6.  We  have  not  proven 
a  convergence  theorem  for  such  scatterers,  but  the  algorithm  seems  to  perform 
quite  well  in  this  case,  and  its  rate  of  convergence  to  be  linear  (see  Table  6.4 
and  Figure  6.5). 

Example  4.1.  Here,  we  reconstruct  a  scatterer  defined  by  the  formula 


q(x)  = 


0.4  if  x  €  [1,2], 
0  otherwise. 


(6.19) 
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Figure  6.4:  Reconstruction  of  Example  2.2  with  a  —  20 


Figure  6.5:  Reconstruction  of  Example  3  with  a  =  10 
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Figure  6.6:  Real  Part  of  the  Scattering  Data  po(&)  in  Example  3 


a 

hk 

Nx 

E2 

t  (sec.) 

10 

0.4 

50 

0.165 

0.230 

20 

0.4 

200 

0.119 

1.51 

40 

0.4 

400 

0.843  xiO-1 

6.03 

Table  6.5:  CPU  Times  and  Accuracies  for  Example  4.1 


In  this  example,  the  scatterer  is  discontinuous,  and  the  conditions  of  Theorems 
5.1,  6.2  axe  violated.  In  fact,  the  integrand  p+  —  p-  does  not  even  converge  to 
zero  as  k  — ►  oo.  The  results  of  this  experiment  are  depicted  in  Figures  6.7,  6.8, 
and  Table  6.5. 


Example  4.2.  In  this  example,  we  reconstruct  a  staircase-shaped  scatterer 
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a 

hk 

Nx 

Ei 

t  (sec.) 

5 

0.2 

o 

o 

t-H 

0.149 

0.430 

10 

0.2 

150 

0.936  xlO'1 

1.18 

20 

0.2 

300 

0.682  xl0_1 

4.40 

Table  6.6:  CPU  Times  and  Accuracies  for  Example  4.2 
defined  by  the  formula 


0 

x  €  (—oo,0.5] 

0.1 

x  €  (0.5, 1.0] 

0.2 

x  e  (l.o,  1.5] 

0.4 

X  €  (1.5, 2.0] 

0.6 

x  €  (2.0, 2.5] 

0.5 

x  €  (2.5, 3.0] 

0.3 

x  €  (3.0, 3.5] 

0.1 

x  €  (3.5, 4.0] 

-0.1 

x6  (4.0, 4.5] 

-0.3 

x  €  (4.5, 5.0] 

-0.2 

x  €  (5.0, 5.5] 

-0.1 

x  €  (5.5, 6.0] 

0 

x  €  (6.0,  oo) 

(6.20) 


This  example  is  similar  to  the  preceding  one,  but  the  shape  of  the  scatterer  is 
more  complicated.  The  results  of  this  experiment  are  shown  in  Table  6.6  and 
Figure  6.9. 

Example  5.  We  investigate  the  sensitivity  of  the  reconstruction  to  perturba¬ 
tions  of  the  initial  data.  In  a  somewhat  crude  test,  we  perturb  the  initial  data  for 
the  algorithm  by  truncating  it  after  a  specified  number  of  decimal  digits  (both 
the  real  and  the  imaginary  parts).  Clearly,  after  such  a  truncation,  the  maxi¬ 
mum  relative  error  is  of  the  order  10D_1  (for  example,  when  the  number  1.999  is 
truncated  after  D  =  1  digits,  the  result  is  1). 

Tables  6.7  and  6.8  demonstrate  the  numerical  results  of  the  reconstruction 
of  Examples  2.1  and  3,  respectively,  with  various  truncations  of  the  input  data. 
In  each  case,  a  was  chosen  sufficiently  large  that  the  error  from  the  truncation 
of  the  trace  formula  due  to  finite  a  (see  (4.194),  (4.198))  is  negligible  compared 
to  the  error  due  to  the  finite  number  D  of  digits  retained.  For  a  given  a,  the 
parameters  hk,  Nx  were  chosen  such  that  accuracy  of  the  reconstruction  was 
not  improved  by  a  further  decrease  of  hk  and/or  increase  of  Nx.  Also  see  Figures 
6.10-  6.14,  comparing  the  scatterers  reconstructed  using  the  perturbed  data  with 
the  prescribed  ones. 
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Figure  6.9:  Reconstruction  of  Example  4.2  with  a  =  20 


D 

a 

hk 

Nx 

E 2 

E°° 

1 

7 

0.1 

100 

0.410 

0.474 

1 

14 

0.1 

300 

0.412 

0.473 

2 

7 

0.1 

100 

0.126 

0.156 

2 

14 

0.1 

300 

0.128 

0.157 

3 

7 

0.1 

100 

0.174  xlO-1 

0.265  xlO'1 

3 

14 

0.1 

300 

0.187  xlO-1 

0.256  xlO"1 

4 

14 

0.1 

300 

0.126  xlO”2 

0.151  xlO'2 

4 

28 

0.05 

600 

0.118  xlO"2 

0.132  xlO'2 

5 

14 

0.1 

300 

0.297  xlO"3 

0.426  xlO'3 

5 

28 

0.05 

600 

0.250  xlO"3 

0.324  xlO-3 

Table  6.7:  CPU  Times  and  Accuracies  for  Example  2.1  with  Truncated  Data 


Figure  6.12:  Reconstruction  of  Example  2.1  with  D  =  3 


D  a  hk  Nx  E 2 


1  10  0.1  150  0.647  0.863 

1  20  0.1  300  0.640  0.852 

2  10  0.1  150  0.121  0.173 


2  20  0.1  300  0.113  0.164 


3  10  0.1  150  0.314  xlO-1  0.602  xlO'1 

3  20  0.1  300  0.206  xlO-1  0.439  xlO'1 


Table  6.8:  CPU  Times  and  Accuracies  for  Example  3  with  Truncated  Data 
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a 

hk 

Nx 

E2 

E°° 

t  (sec.) 

25 

0.2 

250 

0.750  xlO-1 

0.153 

22.6 

25 

0.1 

250 

0.722  xlO"1 

0.145 

45.0 

50 

0.2 

500 

0.512  xlO"1 

0.793  xlO-1 

89.8 

50 

0.1 

500 

0.349  xlO-1 

0.686  xlO"1 

179 

100 

0.1 

1000 

0.158  xlO"1 

0.359  xlO"1 

718 

Table  6.9:  CPU  Times  and  Accuracies  for  Example  6 

We  conclude  the  numerical  examples  by  a  nearly  singular  problem  with  the 
scattering  potential  defined  by  the  formula 

q(x)  =  2  •  e-5(*-O-9*02  -f  sin(5(2x  —  x))  •  sin2(2a:  —  x).  (6.21) 

The  parameters  in  (6.21)  were  chosen  such  that  the  minimum  of  the  function 
1  +  q  is  nearly  zero.  In  fact. 

-0.9953,  (6.22) 

2.10.  (6.23) 

Such  a  scattering  potential  is  extremely  difficult  to  reconstruct  since  the  speed  of 
sound  in  the  scatterer  changes  drastically  (the  ratio  between  the  maximum  and 
minimum  speed  of  sound  is  about  400),  and  the  impedance  functions  p+,p~  have 
large  values,  making  the  ODE  system  (6.3),  (6.4),  (6.5)  stiff.  A  standard  second 
order  Crank-Nicolson  implicit  scheme  was  employed  to  solve  this  problem. 

Example  6.  Table  6.9  demonstrates  the  numerical  results  of  the  reconstruc¬ 
tion  of  the  scattering  potential  (6.21).  Also  see  Figure  6.15  for  the  numerical 
reconstruction. 

The  following  observations  can  be  made  from  Tables  6.1-  6.9  and  Figures 
6.1-  6.15. 

1.  When  the  scatterer  satisfies  the  conditions  of  Theorems  5.1,  6.2,  the  accuracy 
of  the  reconstruction  is  somewhat  better  than  that  predicated  by  these  theorems 
(see  Example  2.1).  This  indicates  that  (as  expected)  the  estimates  (5.13),  (6.12) 
are  somewhat  pessimistic. 

2.  When  the  scatterer  violates  the  conditions  of  Theorems  5.1,  6.2  mildly  (by 
having  discontinuous  derivative  at  the  points  0,  2x),  the  reconstruction  algo¬ 
rithm  still  converges.  Qualitatively,  the  reconstructions  in  Figure  3(a)  should  be 
described  as  good.  A  careful  examination  of  Table  6.4  (and  other  data  not  pre¬ 
sented  in  this  thesis)  shows  that  the  error  of  the  reconstruction  for  such  scatterers 
is  proportional  to  1/a. 


mm,(x)  = 
max,(x)  = 
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Figure  6.15:  Reconstruction  of  Example  6  with  a  =  50 


3.  When  the  scatterer  is  discontinuous  (Examples  4.1,  4.2),  the  algorithm  pro¬ 
duces  results  depicted  in  Figures  6.7,  6.9.  The  oscillatory  behavior  near  the 
discontinuities  resembles  the  well  known  Gibbs  phenomenon.  A  careful  exami¬ 
nation  of  Tables  6.5,  6.6,  (and  other  data  not  presented  in  this  thesis)  shows 
that  in  this  case,  the  point-wise  convergence  is  absent.  In  the  L2-norm,  the  error 
of  the  reconstruction  behaves  like  1  / y/a. 

4.  When  the  initial  data  are  perturbed,  the  resulting  error  of  the  reconstruc¬ 
tion  appears  to  be  proportional  to  the  magnitude  of  the  perturbation,  and  the 
proportionality  coefficient  is  close  to  1.  This  is  a  much  better  estimate  than  the 
one  of  Lemma  5.4  which  bounds  the  condition  number  of  the  algorithm  by  a. 
Qualitatively,  it  can  be  said  that  the  algorithm  is  not  sensitive  to  errors  in  the 
initial  data. 

5.  When  the  speed  of  sound  to  be  reconstructed  varies  by  several  order  of  mag¬ 
nitude,  the  algorithm  encounters  a  mild  difficulty  in  the  form  of  stiffness  of  the 
system  of  equations  (6.3)-  (6.5).  The  problem  is  easily  solved  by  switching  to 
an  implicit  ODE  solver. 
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Chapter  7 


Generalizations  and  Conclusions 


7.1  Generalizations  in  One  Dimension 

Following  is  a  discussion  of  possible  generalizations  of  the  techniques  and  results 
of  this  thesis  in  one  dimension. 

1.  In  their  present  form,  Theorems  5.1,  6.2  require  that  the  scatterer  have  at  least 
four  continuous  derivatives.  Numerical  examples  3-4  of  the  preceding  section 
make  it  abundantly  clear  that  this  is  a  superfluous  requirement.  Obviously, 
Theorems  5.1,  6.2  can  be  generalized  to  at  least  include  the  scatterers  of  the 
type  reconstructed  in  examples  3-4.  Including  the  scatterers  of  examples  3-4 
will  be  somewhat  more  involved,  and  will  require  a  significant  reformulation  of 
Theorems  5.1,  6.2. 

2.  The  algorithm  of  this  thesis  can  be  extended  to  the  Schrodinger  equation. 
The  generalization  is  fairly  straightforward  and  will  be  reported  at  a  later  date. 

3.  In  this  thesis,  we  reconstruct  a  scalar  potential  q  given  the  scattering  data  for 
a  single  Helmholtz  equation.  In  many  problems  of  physical  interest,  the  potential 
has  several  components  (such  as  the  compressional  and  shear  speeds  of  sound  in  a 
medium),  and  the  scattered  data  correspond  to  a  system  of  Helmholtz  equations 
(such  as  equations  of  elastic  scattering,  or  Maxwell’s  equations  in  the  frequency 
domain).  An  extension  of  the  techniques  to  these  cases  appears  to  be  relatively 
straightforward,  and  will  be  reported  at  a  later  date. 

4.  The  impedance  function  formulation  of  the  inverse  problem  (see,  for  example, 
Section  5.2)  can  be  reformulated  as  an  initial  value  problem  for  two  variables 
called  local  reflection  coefficients  R+(x,k )  and  R-(x,k)  defined  below.  Later  in 
Section  7.2,  an  extension  of  this  formulation  to  two  dimensions  will  be  presented. 

For  any  Xo  €  [0,1]  ([0,1]  is  the  support  of  the  scatterer  q),  we  define  the 
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truncated  scatterers  qr(x)  and  qi(x)  by  the  formulas 


9r(x) 


{0  if  x  <  x0, 

q(x)  if  x  >  x0, 


(7.1) 


ii(z)  = 


<?(x)  if  x  <  x0, 
0  if  x  >  xq. 


(7.2) 


Clearly,  the  right-traveling  incoming  wave  <j>inc+(x,k )  =  e'kx  gets  reflected  at 
x  =  xo  due  to  the  discontinuity  of  qr  at  that  point.  Since  the  qT(x)  =  0  for  all 
x  <  xo,  the  total  field,  which  is  the  solution  of  the  Helmholtz  equation,  may  be 
expressed  as  (see  formula  (2.17)  in  Section  2.2) 


<t>+(x,  k )  =  exkx  +  R+(x o,  k)e  ,kx,  for  all  x  <  x0,  (7.3) 


where  R+(xQ,  k)e~tkx  is  the  reflected  (or  back-scattered)  wave.  We  refer  to 
R+(x o,  k )  as  the  local  reflection  coefficient  at  xo- 

Denote  by  p+  and  p+T  the  impedance  functions  for  the  two  scatterers  q  and  qT 
respectively  (see  Section  2.2  for  the  definition  of  the  impedance  functions).  For 
all  x  >  x0,  p+  and  p+r  satisfy  the  same  Riccati  equation  (3.50),  since  qT(x)  = 
q(x).  Furthermore,  they  satisfy  same  initial  condition  (3.52);  consequently 


p+(x,  k)  —  p+ r(x,  fc),  for  all  x  >  x0. 


(7.4) 


Combining  formulas  (2.15),  (7.3),  and  (7.4)  with  the  fact  that  p+r(x,k )  is 
continuous  across  x  =  x0,  we  have  the  formula 


,  ^  1  -  R+{x0,k)e  2ikxo 

P+(x0,  k)  —  1  +  R+{xok)e-2ikx o’ 

(7.5) 

from  which  we  c  'tain 

R+(x0,k)e-^  = 

1  +p+(xo,fc) 

(7.6) 

for  all  real  x0. 

We  observe  that  there  exists  a  real  c  >  0  such  that  for  all  real  x 
k  in  the  upper  half  of  the  complex  fc-plane, 

and  complex 

\R+{x,k)e~2ikx\  <c<l, 

(7.7) 

since  p+  is  uniformly  bounded,  and  Re(p+ )  >  0  is  uniformly  bounded  from  below 
by  a  positive  number  (see  Theorem  4.17  in  Section  4.2). 

It  is  now  clear,  from  the  Riccati  equation  (3.50)  and  formulas  (7.5),  (7.5), 
that  there  is  an  ODE  for  the  local  reflection  coefficient  R+(x,k) 

K+(x,k)  =  ~q(x)  (R+(x,k)e-,tr  + 


(7.8) 
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which  is,  of  course,  a  Riccati  equation  for  R+(x,k)  with  the  variable  k  as  a 
parameter. 

In  a  similar  manner,  the  local  reflection  coefficient  R-(x,k)  is  defined  by  the 
formula 

d>-(x,k)  =  e~tkx  +  R-(xo,k)e'kx ,  for  all  x  >  x0l,  (7.9) 

where  <p~  is  the  solution  of  the  Helmholtz  equation  induced  by  the  left-traveling 


incoming  wave  <^nc_(z,  k)  =  e~tkx  and  the  scatterer  qi  defined  by  formula  (7.2). 
The  following  is  a  list  of  similar  facts  about  f?_. 

For  any  real  x  and  complex  k  in  the  upper  half  plane,  R _  is  connected  with 
p_  via  the  formulas 

,  „  1  -R.(x,k)e2'kx 

P~^X'k)  ~  \  +  R-{X,k)e*kx' 

(7.10) 

R-^ikx  = 

(7.11) 

R _  is  bounded  by  the  formula 

\R~(x,k)e2ikx\  <  c  <  1; 

(7.12) 

and  /?_  satisfies  the  Riccati  equation 

R'_(x,k)  =  jq(X)(R-(x,k)eikx  +  e-'kx)2. 

(7.13) 

Finally,  it  is  easy  to  see  that  there  is  a  trace  formula  associated  with  the  local 

reflection  coefficients  i?+ ,  R_ , 

?'( x )  =  ~^(i  +  ?(*))  (l  +  y/l  +  q{x)j  x 

/OO  /  V 

[R+  (x,  k)e~2tkx  -  R _  (x,  k)e2ikx )  dk.  (7.14) 

-OO  ^  ' 

The  ODE  system  (7.8),  (7.13),  and  (7.14)  for  R+,  i?_,  and  q  is  solved  with 
appropriate  initial  values  of  i?+,  i?_,  and  q  at  x  =  0.  As  is  expected,  the 
convergence  results  are  similar  to  those  presented  in  Section  5.2. 

Numerical  experiments  show  that  the  performance  of  the  algorithm  using  the 
local  reflection  coefficients  is  slightly  and  consistently  better  than  that  using  the 
impedance  functions. 

7.2  Extensions  to  Two  Dimensions 

The  following  is  a  brief  discussion  of  the  generalizations  of  our  inversion  algo¬ 
rithms  in  two  dimensions.  While  the  impedance  function  formulation  of  the 
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inverse  problem  (see  Section  2.2)  can  be  easily  generalized  in  two  dimensions, 
the  reflection  coefficient  formulation  (see  Section  7.1  for  the  1-D  formulation)  in 
two  dimensions  is  not  straightforward.  Below,  we  first  formulate  the  forward  (as 
opposed  to  inverse)  scattering  problem  in  two  dimensions.  We  then  derive  the 
Riccati  equation  for  the  impedance  mappings,  for  the  two  dimensional  inverse 
problem.  A  trace  formula  associated  with  the  impedance  mappings  (similar  to 
the  trace  formula  used  in  Section  5.2)  is  then  presented.  Finally,  for  the  general¬ 
ization  of  the  reflection  coefficient  formulation,  we  will  present  a  Riccati  equation 
for  the  scattering  matrix  in  two  dimensions. 

7.2.1  Forward  Scattering  Problem  in  Two  Dimensions 

We  first  formulate  the  forward  scattering  problem  in  two  dimensions.  Let  us 
consider  the  two  dimensional  Helmholtz  equation 

A <p(x,y)  +  k2  •  (1  +  q(x,y))  ■  <f>(x,y)  =  0.  (7.15) 

where  A  is  the  Laplace  operator.  We  assume  that  the  scatterer  q  :  R2  — >  R  is  a 
smooth  function  and  has  a  compact  support  Li,  and  that  q  >  —  1. 

In  the  forward  scattering  problem,  we  are  interested  in  solutions  of  Helmholtz 
equation  (7.15)  of  the  form 

<P(x,y)  =  Mx,y)  +  iP{x,y),  (7.16) 

where  6  is  referred  to  as  the  total  field,  <po  is  referred  to  as  the  incoming  field, 
and  ip  is  referred  to  as  the  scattered  field.  We  also  say  that  the  total  field  <p  is 
induced  by  the  incoming  field  <f> 0. 

The  incoming  field  is  a  solution  of  the  Helmholtz  equation  (7.15)  with  q  =  0. 
It  can  be  expressed  by  a  linear  combination  of  the  incoming  plane  waves  of  the 
form 

<t>ipw{x,  y)  =  (7.17) 

with  /?  the  direction  in  which  a  plane  wave  travels.  As  can  be  easily  verified,  for 
each  real  /3,  function  (7.15)  is  a  solution  of  the  Helmholtz  equation  (7.15)  with 
q  =  °- 

The  scattered  field  ip  satisfies  the  so-called  Sommerfeld  radiation  condition, 
or  the  outgoing  radiation  condition 

y/r  •  —  i  ■  k  ■  ip^j  — >  0,  cis  r  =  y/x2  +  y2  — >  oo  ,  (7.18) 

and  satisfies  an  inhomogeneous  Helmholtz  equation 

A ip(x,y)  +  k2  ■  (1  +  q{x,y))  ■  ip{x,y)  =  ~k2  ■  q( x,y)  ■  4>0{x,  y).  (7.19) 
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7.2.2  Riccati  Equations  for  the  Impedance  Mappings 

We  define  the  set  of  right-going  plane  waves  by  the  formula  (see  (7.17)  for  the 
definition  of  the  incoming  plane  waves) 


_  f  t/c(x-cos(/3)+y-sin(/?)) 
* rpw  —  \  c 


7T  ^  7T  . 
2  ~  ~  2  s 


(7.20) 


A  solution  of  the  Helmholtz  equation  (7.15)  is  said  to  be  a  right-going  so¬ 
lution,  if  it  is  induced  by  a  right-going  plane  wave,  or  in  general,  by  a  function 
(f>+0  in  the  linear  span  of  the  set  $rpu,  of  right-going  plane  waves  (see  (7.16)). 
Therefore,  the  set  of  right-going  solutions  is  defined  by  the  formula 


$rgs  =  {4>+  =  <t>+ o  +  V>+  I  <P+ o  €  Span($rpuj),t/>  is  the  scattered  field  }.  (7.21) 


We  are  now  ready  to  define  the  impedance  mapping  P+.  For  fixed  x  and  k, 
P+(x,k )  is  a  linear  mapping,  :  Lfoc  —*  Ljoc ,  where  Xfoc  is  the  local  7  2  space  on 
the  real  line.  It  maps,  for  fixed  x  and  k.  a  right-going  solution  <6+  €  $r5s  as  a 
function  of  y  to  its  derivative  with  respect  to  x  as  a  function  of  y,  more  specificly, 

P+(x,  k)  ■  <j>+(x,  y)  =  ^  any  <£+  €  <S>rgs.  (7.22) 

Assuming  that  for  fixed  x  and  k ,  the  right-going  solutions  in  $Tgs  is  dense  in 
Lfoc,  and  combining  (7.22)  with  the  Helmholtz  equation  (7.15),  we  can  easily 
obtain  a  Riccati  equation  for  the  impedance  mapping  P+, 

P;(l,  k)  =  -ik  [Pl(x,  k)  -  -  (/  +  q(x,  •))  ,  (7.23) 

where  I  is  the  identity  operator,  q(x,-)  is  a  diagonal  linear  operator, 

9(^7  •)•  f{y)  =  q(x,y)  •  f{y),  (7.24) 

for  a  fixed  x  and  any  /  €  Lfoc. 

The  impedance  mapping  P_  can  be  defined  by  first  introducing  the  set  of  all 
left-going  solutions  in  a  similar  manner,  which  leads  to  a  Riccati  equation  for  P_, 

PL(x, k )  =  ik(p'(x,k)-±£-2-(I  +  ,(x, •) \  .  (7.25) 

These  Riccati  equations  are  operator  equations;  the  operator  for  example, 
in  the  standard  discretized  form,  is  a  tridiagonal  matrix,  hence  the  term  “matrix 
Riccati  equations”. 
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Combining  the  definition  of  the  impedance  mappings  with  the  standard  WKB 
approximation,  we  obtain  a  trace  formula  (the  derivation  is  omitted  here), 

=  -(1  +  «(*,»))  I"  (P+(i,  k)  -  P.(x,  k))  ■  Fm.[y)ik,  (7.26) 

OX  7T  J—oo 

where  Fone(y)  =  1  for  any  real  y  is  clearly  a  function  in  Lfoc.  Numerical  experi¬ 
ments  show  that  when  the  integral  is  truncated  so  that  it  is  taken  over  the  interval 
[—a,  a],  the  rate  of  convergence  behaves  exactly  like  that  in  one  dimension  (see 
Section  4.3). 

Remark  7.1  While  the  matrix  Riccati  equation  and  the  trace  formula  associ¬ 
ated  with  the  impedance  mappings  are  extremely  similar  to  their  one  dimensional 
counterparts  (see  Section  5.2),  the  spatial  discretization  of  these  two  dimen¬ 
sional  objects  is  quite  different  from  that  in  1-D.  The  fact  that  the  scattered  field 
0  decays  very  slowly  like  l/y/r  makes  spatial  truncation  virtually  impossible. 

In  solving  the  initial  value  problem  of  the  ODEs  (7.23),  (7.25),  and  (7.26), 
the  desired  truncation  in  ^-direction  (see  Remark  7.1)  is  achieved  by  periodizing 
the  scatterer  q  in  y-direction  (the  details  of  this  procedure  are  omitted).  The 
periodized  problem  is  then  truncated  in  frequency  space,  and  discretized  in  a 
procedure  similar  to  that  described  in  Section  6.1.  Numerical  experiments  are 
presently  being  conducted,  and  the  results,  together  with  the  periodizing  proce¬ 
dure,  will  be  reported  at  a  later  date. 


7.2.3  Riccati  Equations  for  the  Scattering  Matrix 

In  a  manner  similar  to  that  in  which  the  impedance  functions  p+  and  p_  are 
connected  to  the  local  reflection  coefficients  R+  and  in  one  dimension  (see 
Section  7.1),  the  two  dimensional  objects  P+  and  P_  are  related  to  linear  map¬ 
pings  called  the  scattering  matrices,  which  are  the  two  dimensional  analogues  of 
the  local  reflection  coefficients  R+  and  R- . 

We  will  define  these  scattering  matrices  in  the  polar  coordinates,  and  will 
present  matrix  Riccati  equations  for  the  scattering  matrices.  The  procedures 
used  here  being  similar  to  those  described  in  Section  7.1,  only  the  main  results 
will  be  presented.  We  will  also  casually  use  several  well  established  results  about 
cylindrical  functions. 

For  any  R  >  0,  we  define  the  truncated  scatterer  qn  by  the  formula 


qFi{r,0)  = 


q{r,  0)  if  r  <  R, 
0  if  r  >  R, 


(7.27) 
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We  first  discuss  the  forward  scattering  from  the  truncated  scatterer  qR.  As  is 
well-known,  for  any  r,  an  incoming  field  (see  Section  7.2.1)  can  be  expressed  as 
the  so-called  Bessel- Fourier  series 

OO 

&ncM)=  £  am- Jm(fcr)e,m®.  (7.28) 

77i=  —  OO 


Is  is  also  well-know  that  for  a  fixed  R,  and  for  any  r  >  R  where  qR(r,0 )  =  0, 
the  scattered  field  induced  by  4>inc  and  qR  can  be  written  as  the  so-called  Hankel- 
Fourier  series 

OO 

</>M)  =  £  /3m(R,k)-Hm(kr)e'm$ ,  (7.29) 

m=— oo 

where  Jm  is  the  first  kind  Bessel  function  of  order  m,  Hm  is  the  first  kind  Hankel 
function  of  order  m.  It  is  well-known  that  once  m  >  | kr,  Jm{kr)  decays  and 
Hm(kr)  grows  like 


Jm(kr^ 

Hm(kr ) 


' e  •  k  ■  r' 


\j2irm  V  2m 


~  —i 
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e-k  ■ 
2m 


(7.30) 

(7.31) 


As  is  well-known,  for  an  incoming  field  of  the  form  (7.28),  there  exists  a 
unique  scattered  field  of  the  form  (7.29),  such  that  <f>  =  <f>tnc  +  tb  is  a  solution  of 
the  Helmholtz  equation  (7.15).  Consequently,  the  linear  mapping 


S(R,k )  ■  {am(R,k),m  =  0,±1,...}  =  {0m(R,  k),m  =  0,±1,...}  (7.32) 


is  well-defined,  and  is  normally  referred  to  as  the  scattering  matrix.  Furthermore, 
the  entries  of  the  matrix  S  decay  very  rapidly,  since  f3m(R,  k)  decays  extremely 
rapidly  due  to  (7.30),  (7.31). 

The  scattering  matrix  S  satisfies  the  Riccati  equation 


Z7TT* 

S'(r,k)  =  —k2-(J(kr)-\-S{r,k)-H(kr))-F-q{r,-)-F~1-(J(kr)  +  H(kr)-S{r,k)), 

(7.33) 

where  J(kr ),  H(kr)  are  diagonal  matrices 


J(kr)  =  diag{J0(kr),  J±i(kr),.. .},  (7.34) 

H(kr)  =  diag{H0(kr),H±i(kr),. . (7.35) 

q(r,  •)  is  a  diagonal  linear  operator 


(7.36) 
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and  F  is  the  Fourier  transform. 

The  derivation  of  this  Riccati  equation  follows  the  procedures  outlined  in 
Section  7.1,  where  Riccati  equations  for  the  local  reflection  coefficients  are  ob¬ 
tained  from  the  Riccati  equations  for  the  impedance  functions.  The  derivation 
is  omitted  here,  and  will  be  reported  at  a  later  date. 

Remark  7.2  While  the  investigation  of  the  inverse  scattering  problem  in  the 
form  of  Riccati  equation  (7.33)  is  in  progress,  numerical  experiments  show  that 
ODE  (7.33)  could  be  useful  for  the  solution  of  the  forward  scattering  problem. 
Unlike  the  Riccati  equations  for  the  impedance  mappings,  where  there  is  a  serious 
problem  with  truncation,  the  scattering  matrix  S  can  be  truncated  easily  since  its 
high-frequency  entries  (those  with  large  indices)  decay  very  rapidly.  We  currently 
use  the  standard  f-th  order  Runge-Kutta  method  to  solve  the  forward  scattering 
problem,  that  is,  starting  from  r  =  0  and  S(r ,  k )  =  0,  solve  the  initial  value 
problem  till  r  =  Rq,  for  some  R0  where  the  entire  circle  contains  the  support  of 
the  scatterer.  The  whole  procedure,  as  is  easy  to  see,  requirs  order  N2 3 4  •  log( N) 
operations  to  obtain  an  N  x  N  scattering  matrix. 


7.3  Conclusions 

An  algorithm  has  been  presented  for  the  solution  of  the  inverse  scattering  prob¬ 
lem  for  the  Helmholtz  equation  in  one  dimension.  The  algorithm  is  based  on  a 
combination  of  the  standard  Riccati  equation  for  the  impedance  function  with  a 
newly  constructed  trace  formula  for  the  derivative  of  the  potential,  and  leads  to 
extremely  accurate  and  efficient  numerical  schemes  for  smooth  scatterers.  The 
principal  differences  between  this  scheme  and  various  layer-stripping  techniques 
(see  [12],  [13],  [14])  are: 

1.  The  algorithm  operates  in  the  frequency  domain,  while  other  efficient  schemes 
are  time-domain  ones. 

2.  While  the  layer- stripping  algorithms  assume  (at  least  conceptually)  that  the 
scatterer  is  piece-wise  constant,  and  are  best  in  this  regime,  our  algorithm  as¬ 
sumes  that  the  scatterer  is  continuously  differentiable.  When  the  scatterer  has  a 
sufficient  number  of  derivatives,  the  algorithm  converges  almost  instantaneously 
(see  Theorems  5.1,  6.2). 

3.  The  principal  drawback  of  the  layer-stripping  algorithms  is  the  fact  that  they 
are  an  essentially  one-dimensional  techniques,  and  the  author  is  not  aware  of 
any  successful  attempts  to  generalize  them  to  higher  dimensions.  Our  techniques 

do  generalize  to  two  and  three  dimensions,  and  in  fact  an  implementation  of  a 
two-dimensional  version  of  the  procedure  is  in  progress. 
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